Jérémie Donà

Jean-Yves Franceschi

Sylvain Lamprier

Patrick Gallinari

PDE-Driven Spatiotemporal Disentanglement

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

PDE-Driven Spatiotemporal Disentanglement

Jérémie Donà, 1 Jean-Yves Franceschi, 1 Sylvain Lamprier, 1 Patrick Gallinari 1,2 1 Sorbonne Université, CNRS, LIP6, F-75005 Paris, France 2 Criteo AI Lab, Paris, France

Motivation

• Disentanglement improves interpretability and prediction • Prior spatiotemporal disentanglement works are often complex and do seldom analysis of its meaning • We aim at grounding spatiotemporal disentanglement on stronger foundations

Separation of Variables as Disentanglement

• Consider the heat equation:

∂u ∂t = c 2 ∂ 2 u ∂x 2 , u(0, t) = u(L, t) = 0, u(x, 0) = f (x), with separable solutions: u(x, t) = µ sin nπ L x φ(x) exp - cnπ L 2 t ψ(t) = ξ φ(x) × ψ(t)
• This method separates factors of variations and can be generalized and applied to numerous systems

Separation of Variables for Observations

• State u v : (x, t) → u v (x, t), observations v = v t 0 , . . . , v t 1 corresponding to a spatial measurement of u: v t = ζ u v (.,
v t = (ζ • U • ξ) φ(•), ψ(t) = D φ, ψ(t)
• In practice, we learn vectorial φ ≡ S and ψ ≡ T t • S and T t 0 are inferred with encoders E S and E T from conditioning frames:

V τ (t 0 ) = v t 0 , . . . , v t 0 +τ

Prediction Constraints

• T ≡ ψ is driven by an ODE: ∂T t ∂t = f (T t ) • Forecasting and alignment losses:

L pred = t v t -v t 2 2 L AE = D S, E T V τ t -v t 2 2

Disentanglement Constraints

• From a strict S invariance constraint to a weaker one to take into account variations of observable content:

∂E S ( V τ (t) ) ∂t = 0 ⇒ L S reg = E S V τ (t 0 ) -E S V τ (t 1 -τ ) 2 2
• Disentanglement loss:

L T reg = T t 0 2 2 = E T V τ (t

  t) • Following the functional separation of variables, we seek φ, ψ, U , ξ and ζ such that: z = ξ φ(x), ψ(t) , u(x, t) = U (z), v t = ζ u(•, t) , with associated ODEs or PDEs on φ, ψ and U • We abstract the unknown spatial coordinates:

0.86 0.7990 0.7713 39.5

										2
										0 )	2
	Input	frames	t=0	Truth	t=6	t=14	t=99	Models	SST Mov. MNIST TaxiBJ Pred. Pred. Swap Pred. MSE SSIM MSE
				MIM				PhyDNet DDPAE	1.27 0.3878 0.2839 41.9 -0.6446 0.6378 -
								MIM	0.91 0.6529 -	42.9
				DDPAE				Ours Ours (no S) 0.95 0.6707 -	43.7
				PhyDNet				Content	input	Input	frames	t=0
				Ours					Swap	Ours	MIM	PhyDNet	Truth	t=4
	Content	input		Swap						t=8

Ninth International Conference on Learning Representations -ICLR 2021, Virtual


Gilles Billen

Josette Garnier

Nitrogen biogeochemistry of water-agro-food systems. The example of the Seine land-to-sea continuum

Keywords: Nitrogen cascade, agro-food system, Seine River, denitrification, eutrophication, leaching, nitrous oxide emission

As an illustration of the usefulness of the concept of water-agro-food systems, this paper describes the mechanisms behind the nitrogen cascade from agricultural soils to the coastal sea rivers through aquifers, riparian wetlands, rivers and streams and the estuary, for the case study of the Seine river watershed and its receiving coastal sea. It is stressed that the structure of the agro-food system, i.e. the way food production and trade, as well as consumption and waste management are organized, is the main determinant of the quality of ground-and surface water as well as of aquatic ecosystem functioning.

Introduction

Nitrogen (N), an essential nutrient in the composition of living organisms, present in proteins and nucleic acids, and also one of the major limiting factors of primary production in terrestrial as well as in aquatic systems, has long been known to play a central role in the control of ecosystem functioning. The massive introduction of reactive nitrogen in the biosphere by human activity since the discovery of the Haber-Bosch process completely changed this control, especially since the 1950s with the intensification and modernization of agriculture to feed the increasing population [START_REF] Smil | Nitrogen in crop production: An account of global flows[END_REF][START_REF] Galloway | Nitrogen cycles: past, present, and future[END_REF][START_REF] Gruber | An Earth-system perspective of the global nitrogen cycle[END_REF]. As early as 1924, Lotka already saw that the development of the nitrogen fertilizer industry represented "nothing less than the ushering in of a new era in the history of the human race, a new cosmic epoch. In the short span of a dozen years -geologically speaking in an instant-man has initiated transformations literally comparable in magnitude with cosmic processes" [START_REF] Lotka | Elements of physical biology[END_REF]. Since then, the multiple environmental and societal threats created by this excessive anthropogenic introduction of reactive nitrogen (Nr) into the biospheric cycle have been underlined, at local, regional and global scales. It affects air quality and greenhouse gas balance, soil quality through acidification, and water quality through nitrate contamination and eutrophication, as well as contributing to biodiversity losses [START_REF] Sutton | The European Nitrogen Assessment: sources, effects and policy perspectives[END_REF]. This paradigm of the nitrogen cascade [START_REF] Galloway | The Nitrogen Cascade[END_REF] can be combined with the vision that close links exist between the biogeochemical functioning of aquatic and terrestrial systems, as well between the diverse aquatic ecosystems succeeding each other along the upstream downstream gradient of complex drainage networks. This land-to-sea continuum view extends the river continuum concept formulated by several authors in the 1980s [START_REF] Vannote | The River Continuum Concept[END_REF][START_REF] Ward | The Four-Dimensional Nature of Lotic Ecosystems[END_REF], but places a special emphasis on the relationship with terrestrial ecosystems, hence on the importance of land cover and land use of the watershed as determinants of the biogeochemistry of rivers and streams including their discontinuities such as reservoirs, estuaries and coastal seas. This view embraces the full land sea continuum, from terrestrial to marine systems.

Another paradigm that was developed during the past decade states that food production can be considered the first determinant of land cover and land use changes, not only through agricultural practices, but also through the structure of the whole system by which food is produced to meet the dietary requirements of the population (see among other studies, the EAT-Lancet Report, Willet et al., 2019). This implies that agro-food systems should be considered as an integral part of hydrosystems for a full understanding of the interactions between human activity and water quality. This is exactly the idea behind the concept of water-agro-food systems as formulated by [START_REF] Garnier | Phosphorus budget in the water-agro-food system at nested scales in two contrasted regions of the world (ASEAN-8 and EU-27)[END_REF][START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]. The purpose of this article is to present the integrated vision we developed of the socioecological metabolism of the Seine watershed territory through the prism of nitrogen transformations and transfers from land to sea, from a long term perspective. We first set the stage by presenting the basic physical geographic features of the Seine basin and the long term trajectory of the agro-food system itself. We then follow the journey of nitrogen from its introduction as a limiting nutrient in cropping systems to its losses to the hydrosphere (via leaching essentially in the form of nitrate, NO 3 -), or to the atmosphere (via emissions of dinitrogen, N 2 ) or nitrous oxide, N 2 O). We also focus on the fate of leached nitrogen through aquifers and riparian zones down to rivers, where it is joined by nitrogen originating from human excretion in urban areas, before being delivered to the sea. The issue of resulting marine eutrophication in the Seine Bay complements the integrative vision. Two contrasting scenarios of what the nitrogen cycle in this land sea continuum could look like in the next 30-50 years are finally presented.

The Seine watershed, its river network and the Seine Bay

The Seine watershed occupies the middle of the large Paris Basin, a sedimentary concentric rock formation, with tertiary limestone at the center, followed by a large cretaceous chalk area separated by a narrow clay band from Jurassic limestones at the periphery. These rock layers, often compared to a pile of plates, lie on an ancient massif of metamorphic and igneous rocks outcropping at the extreme north-east and south-east (Fig. 1a). The Seine River system, with a total watershed area of approximately 75 000 km², drains these geological layers via a regular network of tributaries converging to the large fluvial lower Seine River at the beginning of which Paris lies, 350 km from the outlet to the sea. In the second half of the 19 th century, sluice-locks were installed on the lower Seine and its main tributaries for navigation development. The macrotidal estuarine part of the river begins at Poses, where the last of these dams is located. The estuary has also been largely managed for the purpose of navigation and harbor activities. The Seine estuary delivers its water into the Seine Bight, a rather shallow (25 m mean depth) 3775 km² open bay in the English Channel.

Land use currently comprises 50% cropland, 10% permanent grassland, 0.7% vineyard and orchard, 26% forest and 11% urban and artificialized areas (Fig. 1b). The current population is 17.5 million inhabitants, mostly concentrated in the huge Paris agglomeration and along the lower Seine River sector down to the harbors of Rouen and Le Havre (Fig 1c). The most upstream areas have a declining population.

Long-term trajectory of the agro-food system

In the second half of the 19 th century, like in most regions under temperate climate, the dominant agricultural system in the Seine watershed and in the surrounding regions was based on mixed crop and livestock farming, with crop rotations alternating between one year of forage legume and two years of cereals, and livestock manure being the major fertilizer input in addition to symbiotic fixation of atmospheric nitrogen by legumes (Fig. 2a). This system was able to sustain a rural population of approximately 50 inhabitants per km², based on a diet of 5 kgN/cap/yr with one-third animal proteins, and to still export enough surplus to feed the large Paris agglomeration (Billen et al., 2012a[START_REF] Billen | The Seine Watershed Water-Agro-Food System: Long-Term Trajectories of C, N and P Metabolism[END_REF]data from Le Noë et al., 2018). Gradually during the first half of the 20 th century, and then rapidly during the three decades following World War II, the region in the center of the Paris Basin specialized into stockless crop farming, replacing symbiotic fixation and manure by synthetic fertilizers as the main N inputs to arable soils. Livestock was concentrated in peripheral regions, in the eastern part of the basin (remaining in mixed crop and livestock farming system, Mignolet et al., 2007), and in the "Great-West" (where intensive livestock farming developed, highly dependent on forage import, particularly soybean grain and cakes from South America, Fig. 2b). These trends in territorial specialization of farming systems went hand in hand with a commercial opening of agriculture. Today, the Seine basin area exports 80% of its cereal production, while it imports 70% of the animal proteins consumed by the population (Fig. 2).

The share of animal-to-total protein intake in the human diet, was about 30% at the end of the 19 th century, a proportion consistent with the ratio of animal to vegetal proteins produced by the mixed crop and livestock farming system. Animal protein intake has increased considerably in proportion, reaching 70% of the total protein ingestion in the beginning of the 21st century, and has slowly decreased in recent years (Le [START_REF] Noë | Long-term socioecological trajectories of agro-food systems revealed by N and P flows in French regions from 1852 to 2014[END_REF]. Human excreta in the 19 th century was considered a useful fertilizing resource, collected in cesspools, and then through sewer networks [START_REF] Barles | L'invention des déchets urbains: France 1790-1970[END_REF][START_REF] Esculier | Past and Future Trajectories of Human Excreta Management Systems: Paris in the Nineteenth to Twenty-First Centuries[END_REF]. It is currently treated in wastewater purification plants, which eliminate about 57% of N in excreta as atmospheric dinitrogen through an energetically expensive process of denitrification (using methanol as an electron donor), recover 4% as sludge, only a small part of which is recycled in agriculture, and discharge the remaining part (39%) to surface water [START_REF] Esculier | The biogeochemical imprint of human metabolism in Paris Megacity: a regionalized analysis of a water-agro-food system[END_REF]. The circularity of human excreta management, i.e. the proportion of N in human excreta that is returned to agricultural land, peaked in the 1900's to approximately 40%, and decreased since then to a negligible amount. 1929[START_REF] Constantin | Effects of catch crops, no till and reduced nitrogen fertilization on nitrogen leaching and balance in three long-term experiments[END_REF]-2014(data from Le Noë et al., 2017, 2018).

The consequences of intensification of cropping systems

The trend described above toward territorial specialization of agriculture developed alongside the simplification of crop rotations and intensification of soil N inputs. Crop rotations alternating between forage legumes, wheat and oats, which dominated in the central region of the Paris Basin until the early 1950s, were gradually replaced by rotations alternating between rapeseed, wheat, and barley or grain maize. [START_REF] Lassaletta | 50 year trends in nitrogen use efficiency of world cropping systems: the relationship between yield and nitrogen input to cropland[END_REF] have shown that a robust relationship exists, for a given pedoclimatic context, between yield (Y in kgN/ha/yr, expressed in terms of N content of harvested products over the whole crop rotation cycle) and total N inputs to soil (F in kgN/ha/yr, including manure, symbiotic fixation, atmospheric deposition and synthetic fertilizers over the rotation cycle). This relationship can be represented by a single parameter hyperbolic function of the form: Y = Ymax F/ (F+Ymax), where Ymax characterizes the fertility of the cropping system at maximum N availability, in a given pedo-climatic context. The data gathered from the reconstituted chronicles of crop production and fertilization in the different agricultural regions of the Seine basin over the period 1950-2015 (data from [START_REF] Noë | Long-term socioecological trajectories of agro-food systems revealed by N and P flows in French regions from 1852 to 2014[END_REF] show a trajectory that describes such a relationship, with a relatively constant Ymax, although a trend toward a slightly higher Ymax value (i.e. higher yields at similar or even lower rate of N inputs) is observed after the 1980's (Fig. 3). A similar relationship holds for individual crop rotations in each regions of the Seine basin documented by agronomical surveys for both conventional and organic rotations (Anglade et al., 2015a,b;[START_REF] Rakotovololona | Driving factors of nitrate leaching in arable organic cropping systems in Northern France[END_REF]. Interestingly, the same relationship describes both systems, indicating that there is no intrinsic fertility gap inherent to organic farming practices.

Figure 3. Yield-fertilization relationships for crop rotations in the central regions of the Paris Basin (a)

and in the regions of the eastern periphery of the Seine watershed (b). The ordinate is the annual average export of N imbedded in crop harvest exported over the complete rotation cycle; the abscissa is the annual average N input to the soil as manure, symbiotic fixation, atmospheric deposition and synthetic fertilizer. Gray circles refer to the reconstituted chronicles of the cropping system of Ile-de-France and Picardy (a) and Bourgogne and Lorraine (b) [START_REF] Noë | Long-term socioecological trajectories of agro-food systems revealed by N and P flows in French regions from 1852 to 2014[END_REF]; red and green circles refer to data from current individual crop rotations in conventional (red) and organic (green) farming systems [START_REF] Anglade | Nitrogen soil surface balance of organic vs conventional cash crop farming in the Seine watershed[END_REF][START_REF] Rakotovololona | Driving factors of nitrate leaching in arable organic cropping systems in Northern France[END_REF]. The black curves show the theoretical Y vs F relationship for the stated values of Ymax (in kgN/ha/yr).

Yield-N input relationships are useful, as they allow us to predict the environmental losses of N resulting from agricultural practices. The shape of the relationship implies that these losses defined as the N-surplus, i.e. the difference between N inputs and yieldincreases with intensification, i.e. with increasing fertilization. In other words, nitrogen use efficiency (NUE: the ratio of N yield over N total input) decreases with increasing fertilization. The agricultural N-surplus for cropland of the Seine watershed thus increased from close to 6 kgN/ha/yr in 1850 to as much as 85 kgN/ha/yr at the end of the 1980s, after which it slightly decreased to approximately 50 kgN/ha/yr in 2015 (Fig. 4).

The fate of the surplus, if not stored in the soil organic matter (OM) pool, is either to be lost to the atmosphere or to the hydrosphere. In cropland, leaching generally represents by far the largest fraction of the surplus, unless cover crops are systematically implemented in order to avoid bare soils during the wet autumn and winter season [START_REF] Constantin | Effects of catch crops, no till and reduced nitrogen fertilization on nitrogen leaching and balance in three long-term experiments[END_REF][START_REF] Justes | Réduire les fuites de nitrate au moyen de cultures intermédiaires Conséquences sur les bilans d'eau et d'azote, autres services écosystémiques[END_REF], a practice that also increases soil OM storage [START_REF] Autret | Alternative arable cropping systems: A key to increase soil organic carbon storage? Results from a 16 year field experiment[END_REF][START_REF] Autret | Can alternative cropping systems mitigate nitrogen losses and improve GHG balance? Results from a 19-yr experiment in Northern France[END_REF].

Le Noë et al. (2019a,b) estimated the long term dynamics of organic carbon (C) in the arable soils of France. Sequestration was shown to be mainly the result of soil inputs of humified crop residues, and only secondarily of manure application. In the Seine basin, low rates of sequestration, or even negative rates (i.e. net mineralization) were observed at the beginning of the 20 th century. During the 30 years following World War II, which was marked by a rapid increase in crop productivity, higher rates of C sequestration occurred, which levelled off from the end of the 1980's. As the C:N ratio of soil OM is rather constant at around a value of 10 gC/gN [START_REF] Benoit | A network of organic and conventional crop farms in the Seine Basin (France) for evaluating environmental performance: yield and nitrate leaching[END_REF], Le Noë et al., 2019b), much lower than the C:N ratio of crop residues, the corresponding N sequestration in soil OM can be easily calculated (Fig. 4). As a mean for the arable land of the Seine watershed, it varied from18 kgN/ha/yr in the 1970's to approximately 5 kgN/ha/yr in the 2010s. Another fraction of the surplus is eliminated to the atmosphere through soil denitrification, the process of anaerobic respiration of OM using nitrate as an electron acceptor, occurring in anoxic micro-zones of the soil. Wet conditions, high soil content of nitrate and OM are conditions triggering soil denitrification [START_REF] Tiedje | Denitrification: ecological niches, competition and survival[END_REF][START_REF] Knowles | Denitrification[END_REF][START_REF] Weier | Denitrification and the dinitrogen nitrous-oxide ratio as affected by soil-water, available carbon, and nitrate[END_REF]. Although dinitrogen is the main end product of denitrification, nitrous oxide, a potent greenhouse gas, is also produced [START_REF] Vilain | Nitrous oxide production by nitrification vs. denitrification in experimental conditions[END_REF][START_REF] Benoit | Nitrous oxide production from nitrification and denitrification in agricultural soils: determination of temperature relationships in batch experiments Processes in[END_REF]. For soil denitrification, N 2 O/(N 2 +N 2 O) ratios in the range of 16-40 (median 28) have been estimated (Box 1). Soil N 2 O emission, although of episodic nature, is widely measured using static chambers. [START_REF] Garnier | Longterm changes in greenhouse gas emissions from French agriculture and livestock (1852-2014): from traditional agriculture to conventional intensive systems[END_REF] have proposed an empirical relationship relating soil N 2 O emissions (hence soil denitrification) to climatic data and fertilization rates. Using this relationship, soil N 2 O emission, as well as denitrification, in the Seine watershed was estimated for the period 1900-2015 (Fig. 4). Current values of soil denitrification range at approximately 15-20 kgN/ha/yr.

On an annual basis, leaching is the difference between surplus, N storage in soil OM and denitrification (Billen et al., 2019b). Calculations of the leaching from the surplus minus denitrification and sequestration (Fig. 4), well agrees with on-field measurements, using suction cups, which amounted to 16.2 ± 6.3 kg N ha -1 yr -1 [START_REF] Benoit | A network of organic and conventional crop farms in the Seine Basin (France) for evaluating environmental performance: yield and nitrate leaching[END_REF].

Box 1: N 2 O/(N 2 +N 2 O) ratio of denitrification in soils and aquatic systems.

Considerable disagreement exists in the literature concerning the relative share of N 2 and N 2 O in the final products of denitrification under environmental conditions. Yet, this issue is of considerable practical importance for two reasons. First because N 2 O is a potent greenhouse gas so that the use of denitrification as a way to eliminate nitrate from seepage or wastewater is potentially harmful for the climate, second, because as N 2 O emission is more easy to measure than N 2 production, the former is often used as an indirect way to evaluate the latter.

Many authors report the N 2 O/(N 2 +N 2 O) ratio by denitrification under laboratory conditions where several factors, such as pH, nitrate, oxygen or organic matter concentration are manipulated (e.g. [START_REF] Rochester | Estimating nitrous oxide emissions from flood-irrigated alkaline grey clays[END_REF][START_REF] Hénault | Predicting in situ soil N2O emission using NOE algorithm and soil database[END_REF][START_REF] Billen | Modeling indirect N2O emissions along the N cascade from cropland soils to rivers[END_REF]. The observed effect of these factors may simply reflect the transient perturbation linked to an enzymatic or population readjustment to new conditions imposed by the experimental protocol (Simek et al., 2002a,b). On the other hand, many so called 'in situ' determinations used the acetylene block technique, based on the assumption that acetylene selectively inhibits the last step of denitrification, leaving N 2 O as the sole end product. [START_REF] Butterbach-Bahl | Nitrous oxide emissions from soils: how well do we understand the processes and their controls?[END_REF] have convincingly demonstrated that this assumption is incorrect and that the presence of acetylene modifies other microbial processes, leading to a severe overestimation of the N 2 O/(N 2 +N 2 O) ratio.

Leaving aside the values obtained with these debatable procedures, we are left with only a few direct determinations of the N 2 O/(N 2 +N 2 O) ratio which however show a consistent trend toward one order of magnitude higher ratios in soils (either agricultural or natural) than in wetlands or water bodies (Table B1). The different contribution of nitrification to N 2 O emission beside denitrification can be an explanation for this observation. This also applies to the high N 2 O emission ratio by wastewater treatment facilities. The values in Table B1 have been used in this study to calculate N 2 O emissions from estimated denitrification in river and riparian wetlands and to calculate denitrification from estimated N 2 O emissions from agricultural and forested soils. [START_REF] Kampschreur | Nitrous oxide emission during wastewater treatment[END_REF], only those relating to wastewater treatment involving both nitrification and denitrification were considered.

Nitrogen in the hydrosystem

Nitrate leaching is the main process by which nitrate originating from soil processes is transferred to the hydrosystem. As a first approach, nitrate concentration in the total specific runoff water flux (i.e. precipitation minus evapotranspiration, expressed in mm/yr) over a cropland area can be estimated from the N soil surplus minus denitrification and N sequestration, divided by the total runoff over the same area. A part of the corresponding water flux, which we assume has the same nitrate concentration, rapidly joins surface water (surface or sub-surface rapid runoff), another part infiltrates the vadose zone and the underlying aquifer (aquifer recharge). On a pluri-annual time scale, this partition is termed the base-flow index, the ratio of groundwater (slow) to total (slow and rapid) runoff in the river annual discharge. The presence of large aquifer formations, particularly in the center of the Seine watershed, confers a rather large value to the base-flow index, varying between 0.36 and 0.85 according to the geological substratum (Gallois, pers. comm.), with an average value around 0.70. While rapid runoff reacts rapidly to changes in agricultural practices, base flow responds much more slowly, depending on the time required for groundwater composition to reach equilibrium with the composition of recharge water (Fig. 5). This delay, which varies according to the specific characteristics of the vadose zone and aquifers in the different areas of the watershed, leads to a transient retention of nitrate, which was modelled by [START_REF] Ledoux | Agricultre and groundwater nitrate contamination in the Seine basin. The STICS-MODCOU modeling chain[END_REF] and Beaudoin et al. ( 2016) (Fig. 4). where:

Crforest is the nitrate concentration in leaching water from forest and unmanaged grassland; it has been estimated at approximately 0.5-2 mgN/l from observations of the drainage concentration of natural areas in the Seine basin [START_REF] Billen | River basin nutrient delivery to the coastal sea: assessing its potential to sustain new production of non siliceous algae[END_REF], frcropld is the fraction of land area occupied by cropland, which varies over the long term; Surpl and Acroplnd are, respectively, the calculated N soil surplus of arable land in the considered territory expressed in ktonN/yr, and the area of cropland in ha, as estimated by [START_REF] Noë | Long-term socioecological trajectories of agro-food systems revealed by N and P flows in French regions from 1852 to 2014[END_REF]; Denit is the cropland soil denitrification rate in ktonN/yr estimated from the N 2 O emission calculated by [START_REF] Garnier | Longterm changes in greenhouse gas emissions from French agriculture and livestock (1852-2014): from traditional agriculture to conventional intensive systems[END_REF] taking into account an N 2 O/(N 2 O+N 2 ) ratio of 28 (Box 1); Nsequ is the N sequestration in the soil OM pool of cropland, in ktonN/yr, estimated from the C sequestration rate calculated over the long term by Le [START_REF] Noë | Drivers of long-term carbon dynamics in cropland: a bio-political history (France, 1852-2014)[END_REF], considering a constant C/N ratio of 10 gN/gC (see above); Runoff is the mean long term value of total runoff over the considered area, estimated to be 260 mm/yr, for the whole watershed (63000 km² at Poses), based on the long term record of discharge of the Seine river at the outlet (Poses, 512 m 3 /s in average over the last 50 yrs).

It can be seen that the delay in the response of aquifers to the increase in nitrate concentration of recharge water is typically of approximately 10-25 years (Fig. 5). This is in good agreement with determinations of groundwater age using measurements of dissolved chloro-fluoro-carbon compounds [START_REF] Vergnaud-Ayraud | La datation des eaux souterraines par analyse des CFC : un outil de gestion durable de la ressource en eau[END_REF] in several drinking water wells in the Seine basin [START_REF] Anglade | Reconquérir la qualité de l'eau en régions de grande culture : agriculture biologique et reconnexion avec l'élevage[END_REF]Beaudoin et al., pers. comm.).

Figure 5. Long term trends in groundwater nitrate concentrations from data on several individual

springs in the central Seine Basin (Vulaines-lez-Provins and Saint Brice, Aquibrie, pers. comm.;Melarchez, Anglade et al., 2017) with varying residence time, compared with the calculated concentration of recharge water in Ile-de-France, calculated by relation [START_REF] Esculier | The biogeochemical imprint of human metabolism in Paris Megacity: a regionalized analysis of a water-agro-food system[END_REF].

Both base flow and surface runoff, with their own differing composition, have to cross the riparian zone before reaching surface water, unless it has been by-passed by the installation of tile drains, a widespread practice in some regions of the Seine basin until the 1980s and is still the case for almost 11% of the Seine watershed area. The riparian zone generally consists of a more-or-less wide area of uncultivated wetlands, where ground-and subsurface waters flowing from the watershed reach the upper, biogeochemically active layer of the soil (Fig. 6). The soil denitrification potential is expressed there by its capacity to reduce the nitrate loading of the water flow, before it reaches the river itself. Based on this conceptual scheme and on a basin-wide delimitation of riparian zones resting on land use and on topographic indexes [START_REF] Berthier | Enveloppes des milieux potentiellement humides de la France me´tropolitaine. Notice d'accompagnement. Programme de modélisation des milieux potentiellement humides de France[END_REF], Billen et al. (2018b[START_REF] Billen | Modeling indirect N2O emissions along the N cascade from cropland soils to rivers[END_REF] modeled nitrate and N 2 O dynamics associated with these riparian processes (Fig. 6), and provided an estimate of their significance at the local and basin scales (Fig. 4). These data show that although riparian areas are hotspots of denitrification and N 2 O emission in the sense that per area flux values are generally higher there than at other sites in the watershed, the integrated value at the watershed scale is lower, because of the small proportion of watershed area covered by active riparian wetlands (see [START_REF] Garnier | Curative vs. preventive management of nitrogen transfers in rural areas: lessons from the case of the Orgeval watershed (Seine River basin, France)[END_REF]. Riparian denitrification currently amounts about 15% of the diffuse flux of nitrogen at the watershed scale. This figure could have been twice higher in earlier times, before riparian wetlands were drained and/or converted to cropland. The same holds true for in-stream denitrification processes which, although significant when expressed in terms of local per area values, affect nitrogen dynamics only marginally compared with soil and riparian processes. However, downstream from the largest point of wastewater discharge, microbial activity occurs, either heterotrophic mineralization of organic N into ammonium, or nitrification of ammonium, and locally affects water quality and N 2 O emissions [START_REF] Garnier | Nitrous oxide (N2O) in the Seine river and basin: observations and budgets[END_REF]. In this respect, the sector of the Seine River downstream from the discharge of the Seine-Aval wastewater purification plant treating the largest part of the Paris effluents, plays a significant role. The never-ending race between an increasing N load collected by the sewers to the wastewater plants, and the improvement of the treatments applied is reflected in the long-term series of water quality observations and modelling results of water quality along longitudinal profiles of the lower Seine River course, as shown by [START_REF] Aissa-Grouz | The response of river nitrification to changes in wastewater treatment (The case of the lower Seine River downstream from Paris)[END_REF], [START_REF] Romero | Long-term water quality in the lower Seine: lessons learned over 4 decades of monitoring[END_REF] and Garnier et al. (2019c).

Denitrification, and N 2 O emissions, is also significant in the large reservoirs fed by the upstream course of the Seine river and its tributaries the Aube and Marne river that were established in for the purpose of discharge regulation. During the summer period of water stagnation in these reservoirs as well as in farm ponds, nitrate concentration is strongly reduced before restitution to the river network [START_REF] Garnier | Ecological functioning of a large reservoir in the upstream basin of the river Seine (Marne reservoir, France)[END_REF]2014;[START_REF] Passy | Restoration of ponds in rural landscapes: Modelling the effect on nitrate contamination of surface water (the Seine River Basin, France)[END_REF].

N delivery to the coastal zone and marine eutrophication

The flux of N delivered by the Seine River at the outlet of the river network (more specifically at Poses, the entrance of the tidal estuarine sector) reflects all the diffuse and point sources distributed in the watershed area, as well as the various processes leading to the transformation, immobilization, and transfer of N loading along the land-water continuum [START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]. Figure 4 summarizes the long term trends of these sources and sinks and shows how the flux in N delivery has increased as the result of these trends.

Overall, the latter currently amounts to 42-49% of net N sources to the watershed. N "retention" consists of 8% soil storage, 15% soil denitrification, 10% riparian denitrification Because of the heavy morphological transformations undergone by the Seine estuary for the purpose of facilitating navigation, the N retention function of this system has been considerably reduced. The flux of N delivered to the estuary is now transferred to the Seine Bay, with only 13% (about 20 ktN/yr) being retained on an annual basis by denitrification [START_REF] Romero | Modeling the biogeochemical functioning of the Seine estuary and its coastal zone: Export, retention, and transformations[END_REF].

In the Seine Bay on the other hand, the input of nitrogen, together with the other nutrients phosphorus and silica, determine the state of eutrophication of the ecosystem. [START_REF] Billen | River basin nutrient delivery to the coastal sea: assessing its potential to sustain new production of non siliceous algae[END_REF] defined the N-ICEP, an indicator of the potential role of nutrient riverine delivery in promoting the growth of undesirable, non-siliceous marine phytoplankton species in coastal bights. The concept behind this indicator is that when nutrients are brought into the sea in proportions corresponding to the Redfield ratio of algal composition, they promote the development of diatom communities able to feed a linear food chain leading to fish, while if N and/or phosphorus are in excess over silica with respect to diatoms requirements, algal blooms dominated by potentially toxic or harmful species develop, including dinoflagellates or even low-silicified diatoms such as Pseudo-nitzschia, producing domoïc acid in response to silica depletion [START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]. ) is able to reasonably simulate diatom and nanoflagellate blooms [START_REF] Passy | Eutrophication modelling chain for improved management strategies to prevent algal blooms in the Seine Bight[END_REF][START_REF] Romero | Modeling the biogeochemical functioning of the Seine estuary and its coastal zone: Export, retention, and transformations[END_REF][START_REF] Desmit | Reducing marine eutrophication may require a paradigmatic change[END_REF], including Pseudo-nitzschia [START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF], in response to the reconstructed riverine nutrient fluxes. [START_REF] Billen | River basin nutrient delivery to the coastal sea: assessing its potential to sustain new production of non siliceous algae[END_REF]; open circles: values calculated with the Seneque-Riverstrahler model for actual conditions from 2002 and 2014 [START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]

Two contrasting scenarios for the future

To a large extent, the current knowledge on nutrient cycling in the past and present wateragro-food system also allows us to predict the environmental effects of possible future changes in the design of human food provision. For example, two contrasting scenarios of the future of French agro-food systems have been constructed (Billen et al., 2018a;[START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]. One of these scenarios, called "Opening and Specialization" (O/S) is constructed by extending the historical trends observed in the French agro-food system during the past few decades: livestock and permanent grasslands disappear from the central regions of the watershed but livestock farming intensifies in the peripheral regions with an increased import of soybean; cereal production for export keeps growing; a human diet rich in animal protein persists; current environmental regulations, including those related to mineral and organic N fertilization, are maintained. The second scenario, called "Autonomy, Reconnection and Demitarian Diet" (A/R/D), assumes the generalization of organic farming practices with long and diversified crop rotations and no industrial fertilizers, local livestock and crop farming reconnection, and a human diet of half the current content of meat and milk. At the scale of France, both scenarios were shown to be able to meet the needs of the projected French population for 2050, in both cases with ample exportation of cereals. In the O/S scenario France also exports animal products but imports massive amounts of feed; by contrast in the A/R/D scenario no significant net export of meat and milk occurs and no import of feed is required. At the scale of the Seine basin (which hosts 25% of the French population), the A/R/D and O/S scenarios require importing of, respectively, 33% and 86% of animal protein requirements (Fig. 8a, Table 2).

Table 1 : Main functional characteristics and production fluxes of the agro-food system of the Seine watershed (in bold) and the whole of France (in italic) for the present situation (2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014) and for the two prospective scenarios O/S and A/R/D (Billen et al., 2018a). (Billen et al., 2018a). b. Simulated nitrate concentration in the Seine drainage network in both scenarios compared with the present situation (calculated with the Seneque-Riverstrahler model, [START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]. c. calculated N 2 O emission by agricultural soils in both scenarios compared with the present situation [START_REF] Garnier | Longterm changes in greenhouse gas emissions from French agriculture and livestock (1852-2014): from traditional agriculture to conventional intensive systems[END_REF] These works underline that radical changes in the structure of the agro-food system are required to reach a good ecological status of the hydrosystem, rather than simple improvements in farming practices or policies to increase N use efficiency during the different transformation stages. Climate change is not taken into account in these scenarios which typically cover the next 30 years.

2004-2014

ARD scenario OS scenario Reference

Conclusion

During the past 30 years, the field of biogeochemistry, which was still in its infancy at the time the eponymous journal was founded, has evolved toward a major discipline at the core of several issues of vital importance for human life on Earth. Numerous paradigms have been formulated during this period, which today serve as the backbone of environmental sciences. The river continuum concept was extended to a full land-to-sea continuum vision; combined with the synoptic approach of agro-food systems it resulted in the water-agro-food concept that we applied here to the study of the nitrogen cascade.

The case study of the Seine basin illustrates how the knowledge gained during the past 30 years of the long term PIREN-Seine program (since 1989), on the processes involved in the N cycle and its perturbations, has made it possible to quantitatively describe the associated fluxes at the scale of a whole watershed, over a timespan of more than one century.

Changes in land use of the watershed, and particularly in the agricultural activities and in the management of urban wastewater resulting from human nutrition are by far the main driver of the observed trajectory of the functioning of the water-agro-food system. Several possible future scenarios for the agro-food system can be constructed and their consequences estimated in terms of water quality and ecological functioning of the whole continuum of aquatic systems from rivers to the sea; their impact on the climate system through greenhouse gas emission can also be calculated. What is more challenging is to incorporate the effect of global climate change into such prospective scenarios at the regional scale. If robust estimations of the future changes in temperature, rainfall, and runoff are available for the next 30-50 years, it is feasible to incorporate these changes and calculate how they affect the control on the kinetics of the processes described in the ecosystem models. Such explorations have already been initiated with changes in hydrology but not yet with temperature [START_REF] Raimonet | Coastal eutrophication potential under future climate change: a landward perspective[END_REF][START_REF] Garnier | Nutrient inputs and hydrology together determine biogeochemical status of the Loire River (France): current situation and possible future scenarios[END_REF]. A still more challenging issue is however to take into account the indirect effects of climate change on agriculture. Current cropping systems may well be severely affected by these changes, and new crops and farming practices might be introduced with consequences that are difficult to predict using our current knowledge. This issue probably remains a research priority for the coming decades.

Figure 1 .

 1 Figure 1. a. The Seine River system and the geological formations of its watershed and surrounding areas (source: BRGM www.brgm.fr). From left to right b. Current distribution of land use (source: Corine Land Cover 2012, www.data.gouv.fr/fr/datasets/corine-land-cover-occupation-des-sols-enfrance) and land cover variations in the Seine watershed since 1850 (data from Le Noë et al. subm) c. Current distribution of population density (INSEE 2017 https://www.insee.fr/fr/information/2008354.) and population evolution in the Seine watershed (data from Le[START_REF] Noë | Long-term socioecological trajectories of agro-food systems revealed by N and P flows in French regions from 1852 to 2014[END_REF]

Figure 2 .

 2 Figure 2. GRAFS representation of the agro-food system of the Seine watershed in the periods of

Figure 4 .

 4 Figure 4. Nitrogen budget of the Seine river watershed(limited at Poses, 69000 km²). Positive fluxes are the inputs to the water system :[START_REF] Esculier | The biogeochemical imprint of human metabolism in Paris Megacity: a regionalized analysis of a water-agro-food system[END_REF] cropland N surplus defined as the share of N soil inputs that is not exported by crop harvest; (2) N inputs to the hydrosystem from non-cropland soils (3) urban wastewater discharge to the river system after treatment. Negative fluxes refer to "retention" processes:[START_REF] Esculier | The biogeochemical imprint of human metabolism in Paris Megacity: a regionalized analysis of a water-agro-food system[END_REF] N sequestration in the soil organic matter pool, (2) soil denitrification; (3) temporary retention in vadose zone and aquifer; (4) riparian denitrification and (5) in-stream denitrification. The black curve refers to the river N delivery at the outlet of the river network: it roughly corresponds to the difference between inputs and retentions.

  Figure 5 illustrates this response by comparing the reconstituted time variations of nitrate concentration in the recharge water with the observed concentration in a number of well-documented springs in the central part of the Paris Basin. The mean nitrate concentration in the recharge water (CrNO3, mgN/l) of the central Paris Basin (Ile-de-France region) was calculated for the period 1850-2015 from the following relationship: CrNO3 = Crforest * (1-frcropld) + [(Surpl-Denit-Nsequ)*10 9 / Acropld ] / (Runoff*10)] * frcropld (1)

Figure 6 .

 6 Figure 6. Block diagram showing the transfer of nitrate from watershed soils to river channel through riparian wetlands, and the associated denitrification and N 2 O emission. In-stream and reservoir denitrification processes are also represented.

  The long-term trend of the N-ICEP in the flow delivered to the Seine Bight is shown in Figure 7. While negative N-ICEP values would occur under pristine conditions, positive values have been calculated since 1850; they increased rapidly after 1950 and are still positive today, although a slow decreasing trend is observed since the 1980s. The ECOMARS-3D model of the Channel developed by IFREMER (http://wwz.ifremer.fr/mars3d/Le-modele/Descriptif/, Menesquen et al., 2018

Figure 7 .

 7 Figure 7. Long-term trends of N-ICEP in the nutrient flow delivered by the Seine River into the Seine estuarine sector (at Poses). Black circles and bars: values calculated for median, wet or dry hydrology[START_REF] Billen | River basin nutrient delivery to the coastal sea: assessing its potential to sustain new production of non siliceous algae[END_REF]; open circles: values calculated with the Seneque-Riverstrahler model for actual conditions from 2002 and 2014[START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]

  Simulated water quality in the Seine drainage network shows much lower nitrate contamination in the A/R/D than in the O/S scenario (Fig.8b), and a much lower N loading to the Seine Bight (respectively 3200 and 2360 kN/km²/yr for A/R/D and O/S, vs. 4200 kN/km 2 / yr for 2010-2014). N 2 O emission from cropland soils is also considerably lower in the A/R/D than in the O/S scenario (Fig.8c).

Figure 8

 8 Figure 8. a. GRAFS representation of the agro-food system of the Seine watershed in the O/S and A/R/D scenarios for the year 2050(Billen et al., 2018a). b. Simulated nitrate concentration in the Seine drainage network in both scenarios compared with the present situation (calculated with the Seneque-Riverstrahler model,[START_REF] Noë | Carbon dioxide emission and soil sequestration for the French agro-food system: present and prospective scenarios[END_REF]. c. calculated N 2 O emission by agricultural soils in both scenarios compared with the present situation[START_REF] Garnier | Longterm changes in greenhouse gas emissions from French agriculture and livestock (1852-2014): from traditional agriculture to conventional intensive systems[END_REF]

Table B1 :

 B1 In situ determination of the N 2 O/(N 2 +N 2 O) ratio in soils, wetlands and water bodies, excluding those obtained using the acetylene block technique.

	Environment	N 2 O/(N 2 +N 2 O) ratio, % nb measurements	reference	.
	Agricultural soils		28 ± 12	30	Schlesinger, 2009
					Butterbach-Bahl et al, 2013 *
	Natural soils		28 ± 17	18	ibid
	Wetlands, flooded soils		4 ± 5	14	ibid
	Rivers		0.9 ± 0.4	40	Beaulieu et al., 2011
	Rivers and riparian zones	1.5 ± 0.5	3	Billen et al, 2020
	Wastewater trtmt facilities**	17 ± 10	15	Kampschreur et al., 2009
	* The data compiled by Schlesinger, 2009 (Suppl Mat) and Butterbach-Bahl et al 2013 (Suppl Mat) have been
	combined after deleting those in the former publication that were obtained through the acetylene block technique.
	** From the data compiled by

  and 16% in stream denitrification. The largest part of 'retention' processes thus is attributable to atmospheric emission, partly under the form of N 2 O contributing to climate change.

	Watershed soils
	N 2 , N 2 O	N 2 , N 2 O
			Riparian
			wetlands
	Aquifer
		N 2 , N 2 O	N 2 , N 2 O
	River
			Reservoirs

Acknowledgments

This work was carried out within the scope of the PIREN-Seine program, supported by the Seine-Normandie Water Agency and several other partners. Several others projects are also acknowledged, especially, the NEREIS and RESET projects (GIP-Seine-Aval), as well as FLAM (Liteau & Seine-Normandie Water Agency), ESCAPADE (ANR), HYDROGES (ADEME) and ABAC (Région Ile-de-France, Seine-Normandie Water Agency and Eau de Paris). It also benefited from the scientific context of the FIRE (FR 3020), a federation for environmental research in the Ile-de-France Region and of the eLTER Zone Atelier Seine, promoted by the French National Center for Scientific Research (CNRS-INEE).


Beniamin Bogosel
email: beniamin.bogosel@cmap.polytechnique.fr

Regularity result for a shape optimization problem under perimeter constraint

Keywords: eigenvalues, shape optimization, regularity AMS Subject Classifications. 49Q10, 35J25

We study the problem of optimizing the eigenvalues of the Dirichlet Laplace operator under perimeter constraint. We prove that optimal sets are analytic outside a closed singular set of dimension at most d-8 by writing a general optimality condition in the case the optimal eigenvalue is multiple. As a consequence we find that the optimal k-th eigenvalue is strictly smaller than the optimal (k + 1)-th eigenvalue. We also provide an elliptic regularity result for sets with positive and bounded weak curvature.

Introduction

For a given domain Ω ∈ R d we can consider the eigenvalue problem for the Laplace operator with Dirichlet boundary conditions defined by

-∆u = λ k (Ω)u in Ω u = 0 on ∂Ω.
This operator has compact resolvent if L 2 (Ω) injects compactly in H 1 (Ω) and in this case the spectrum consists of an increasing sequence of eigenvalues

0 < λ 1 (Ω) ≤ λ 2 (Ω) ≤ ... ≤ λ k (Ω) ≤ ... → ∞.
The optimization of the eigenvalues of the Dirichlet Laplace operator under various constraints is a classical problem which recently had important developments. Authors were initially concerned with considering volume constraints, but in some recent works like [START_REF] Bucur | Minimization of λ 2 (Ω) with a perimeter constraint[END_REF], [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] the perimeter constraint gained some interest. In fact, having a perimeter constraint allows the use of techniques related to perimeter quasi-minimizers [START_REF] Tamanini | Variational problems of least area type with constraints[END_REF] which allow important gains in regularity properties. More precisely, the basic form of the problem which we consider is the minimization of λ k (Ω) under a perimeter constraint min

Per(Ω)=c λ k (Ω).
Using the nice behavior of these eigenvalues under homotheties, λ k (tΩ) = λ k (Ω)/t 2 one can deduce that the above problem is equivalent, up to a homothety, to the problem min

Ω∈R d (λ k (Ω) + Per(Ω)) . (1)
This latter formulation has the advantage of being unconstrained and will be used throughout the paper. In [START_REF] De Philippis | Regularity of minimizers of shape optimization problems involving perimeter[END_REF] the authors study existence and regularity properties of problems of the form where G is the Dirichlet energy or a spectral functional. The authors focus on showing how to recover the C 1,α regularity of the relaxed optimizers. They also describe how to perform a bootstrap argument in certain cases and be able to show that optimal sets are C ∞ . The case where the functional depends on spectral quantities is not treated in full detail as one needs to take care of the case where the optimal eigenvalues are multiple, thus losing some differentiability properties.

The aim of this paper is to improve the known regularity properties for problem [START_REF] Ambrosio | Functions of bounded variation and free discontinuity problems[END_REF]. In fact for the same problem with volume constraint regularity issues are more delicate. Existence results are proved in [START_REF] Bucur | Minimization of the k-th eigenvalue of the Dirichlet Laplacian[END_REF], [START_REF] Mazzoleni | Existence of minimizers for spectral problems[END_REF] in the class of quasi-open sets, which are basically level sets of H 1 functions. Initial regularity results are presented in [START_REF] Bucur | Lipschitz regularity of the eigenfunctions on optimal domains[END_REF] where the authors prove that in certain cases the solution to the analogue spectral optimization problem under volume constraint is an open set. From now on we consider only the perimeter constraint.

The Faber-Krahn inequality and the isoperimetric inequality immediately imply that the first eigenvalue is minimized by the ball in any dimension. A first result concerning higher eigenvalues is due to Bucur, Buttazzo and Henrot [START_REF] Bucur | Minimization of λ 2 (Ω) with a perimeter constraint[END_REF]. The authors consider the minimization of the second eigenvalue among domains Ω in R 2 . Some results particular to this case, like convexity of the optimizers, H 2 regularity for elliptic problems on convex domains [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF] and simplicity of the second eigenvalue, allow the use of a bootstrap procedure which shows that the optimal set is C ∞ .

In [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] De Philippis and Velichkov show that problem [START_REF] Ambrosio | Functions of bounded variation and free discontinuity problems[END_REF] has solutions in every dimension and the optimal shapes have regularity at least C 1,α outside a closed set of Hausdorff dimension less or equal than d -8. The fact that minimizers are not known and do not have a simple structure even in the case of the second eigenvalue [START_REF] Bucur | Minimization of λ 2 (Ω) with a perimeter constraint[END_REF] motivated the numerical study of problem [START_REF] Ambrosio | Functions of bounded variation and free discontinuity problems[END_REF]. Simulations were performed in R 2 and R 3 with two distinct methods in [START_REF] Antunes | Optimisation of eigenvalues of the Dirichlet Laplacian with a surface area restriction[END_REF] and [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF]. In addition to numerical simulations Bogosel and Oudet [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF] also present an optimality condition applicable in the case of multiple eigenvalues inspired by methods from [START_REF] Soufi | Domain deformations and eigenvalues of the Dirichlet Laplacian in a Riemannian manifold[END_REF]. The need of these more general optimality conditions is motivated by the numerical observations. Indeed, we observe that except some particular cases the eigenvalues of the optimal shapes seem to have multiplicity higher than one. Thus classical methods based on optimality conditions obtained directly from the derivative of the eigenvalue, like the one used in [START_REF] Bucur | Minimization of λ 2 (Ω) with a perimeter constraint[END_REF], are not always applicable. The more general optimality condition presented in [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF] states that if the optimal set Ω * is of class C 3 then there exist a family of eigenfunctions {φ i } m i=1 such that the mean curvature H of Ω * can be expressed as

H = m i=1 (∂ n φ i ) 2 .
This optimality condition allows the use of a bootstrap argument which shows that optimal shapes are smooth. These optimality conditions were also used in [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF] as a tool for validating the numerical computations.

The purpose of this note is to remove the additional C 3 regularity hypothesis used in the deduction of the above optimality condition. We prove the following result. The arguments used in the proof are different from the ones used in [START_REF] Bourdin | Optimal partitions for eigenvalues[END_REF] which used heavily the C 3 regularity. There are also common points such as the analytic parametrization of eigenvalues under regular perturbations and the use of a variant of the Hahn-Banach theorem to conclude. A summary of the steps employed in the proof is presented below.

1. We want to be able to recover some differentiability information in the case the eigenvalue is multiple. It is well known (see [START_REF] Henrot | Variation et optimisation de formes[END_REF]) that in this case the eigenvalue may not be differentiable. Nevertheless, applying a classical result of Rellich [START_REF] Rellich | Störungstheorie der Spektralzerlegung[END_REF] we can deduce that the eigenvalues and eigenfunctions of an analytic perturbation of the domain can be reparametrized analytically. This fact implies, in particular, that even if we lose differentiability we still recover some facts about directional derivatives.

2.

A second aspect is writing the lateral derivatives of the eigenvalues using expressions similar to the shape derivative formula. In order to do this we would like to have more information about the regularity of the eigenfunctions. In the work of Grisvard [14, Theorem 3.

2.1.2] it is proved that if Ω is convex then solving -∆u = f in H 1 0 (Ω) with f ∈ L 2 (Ω) implies that u is a H 2 (Ω) function.
The proof of this uses the fact that convex sets can be approximated uniformly from inside with C 2 convex sets. This brings us naturally to the next point. [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] it is proved that solutions of problem (1) exist, are of class C 1,α and that any optimizer has non-negative and bounded weak mean curvature. The question which arises is the following: can we approximate C 1,α sets with weak non-negative bounded mean curvature with smoother sets (at least C 2 ) with non-negative mean curvature? The answer is affirmative and a proof is presented in [START_REF] Metzger | No mass drop for mean curvature flow of mean convex hypersurfaces[END_REF]Lemma 3.8] 4. The final step consists in using a variant of the Hahn-Banach theorem which allows us to deduce that the mean curvature is a convex combination of squares of normal derivatives of an orthonormal basis of the eigenspace corresponding to the k-th eigenvalue. Once we have this formula we use a result of Landais [START_REF] Landais | A regularity result in a shape optimization problem with perimeter[END_REF] on the regularity of solutions to the mean curvature equation to deduce that the optimal set is C 2,α . Results concerning the regularity of solutions to equations of the mean curvature type can also be found in [START_REF] Ambrosio | Functions of bounded variation and free discontinuity problems[END_REF]Section 7.7]. A bootstrap algorithm similar to the one used in [START_REF] Bucur | Minimization of λ 2 (Ω) with a perimeter constraint[END_REF] or [START_REF] Landais | A regularity result in a shape optimization problem with perimeter[END_REF] allows us to conclude that optimal shapes are C ∞ .

In

After proving the regularity result we turn our attention to the case where the eigenvalues of the optimal set are multiple. We prove a result which was observed in the numerical computations from [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF]. The result says that when the optimal k-th eigenvalue is multiple then the multiplicity cluster ends at λ k , namely ... = λ k-1 (Ω * ) = λ k (Ω * ) < λ k+1 (Ω * ). This implies, in particular, that the optimal costs

c k = min Ω∈R n (λ k (Ω) + Per(Ω)) , d k = min Per(Ω)=1 λ k (Ω) (2)
form a strictly increasing sequence.

Elliptic regularity for sets of positive mean curvature

As the results of [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] prove that the solutions of problem (1), denoted by Ω * in the sequel, are C 1,α domains (outside a residual closed set of dimension at most d -8) with non-negative weak mean curvature, we start from this setting and we prove the H 2 regularity of the eigenfunctions for the optimizers of (1). In the first part of this section we restrict ourselves to dimension d at most 7, where we do not have singularities. In this case we conclude that u ∈ H 2 (Ω). At the end of the section we give an argument which allows us to deduce local H 2 regularity outside the singular set. This type of results are new up to the author's knowledge. Results relating the regularity of the solution to the curvature of the boundary can be found, for example in [START_REF] Gilbarg | Elliptic partial differential equations of second order[END_REF]Section 14.3], however the problem we are interested in is not covered by these. The results in [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF]Chapter 3] concerning elliptic problem on convex domains are similar to the ones we are interested in and we follow similar methods in the proofs below.

Thus, following the results of [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] we know that solutions of (1) can be represented locally around each point in the regular part of the boundary as epigraphs of C 1,α functions. Thus, in a local coordinate system with origin in a regular point x 0 , the boundary ∂Ω * can be represented as the graph of h : D ⊂ R d-1 ∈ R, where D is an open set containing 0. Moreover, following the results of [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF], for every ψ ∈ C 1 (D), ψ ≥ 0 we have 0

≤ D ∇h • ∇ψ 1 + |∇h| 2 ≤ K ∞
where K ∞ is a constant. Thus the mean curvature of an optimizer Ω * can be represeted as an L 2 (∂Ω) integrable function, since Ω * has finite perimeter. Note that if Ω * is C 2 then these inequalities simply say that its mean curvature is non-negative and bounded above by K ∞ .

In [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF]Chapter 3] we can find results concerning the regularity of elliptic problems in terms of the mean curvature of the boundary. In particular, for convex sets the solution of -∆u = f in H 1 0 for f ∈ L 2 is an H 2 function. The scope of this section is to extend such results to the case of C 1,α domains which have weak non-negative mean curvature. In order to do this we follow a similar approach. In a first stage we approximate C 1,α domains with non-negative weak mean curvature with C 2 sets with non-negative mean curvature. This approximation stage is presented in [START_REF] Metzger | No mass drop for mean curvature flow of mean convex hypersurfaces[END_REF] and uses notions related to weak mean curvature flows. A second stage, which is a simple consequence of Theorem 3.1.2.1 and Remark 3.1.2.2 in [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF] gives some results for C 2 domains with non-negative mean curvature. We conclude the general case using an approach similar to Theorem 3.1.2.1 from [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF].

Proposition 2.1. Let Ω be a C 1,α domain with non-negative weak mean curvature in L 2 (∂Ω). Then for each ε > 0 there exists a smooth set Ω ε , with non-negative mean curvature such that

Ω ε ⊂ Ω and d H (∂Ω, ∂Ω ε ) < ε (where d H represents the Hausdorff distance).
Proof: We see that the set Ω verifies the hypotheses of Lemma 3.8 in [START_REF] Metzger | No mass drop for mean curvature flow of mean convex hypersurfaces[END_REF]. Indeed, Ω has nonnegative weak mean curvature in L 2 (∂Ω) so there exists a sequence (Ω ε ) ε>0 such that Ω ε → Ω in C 1 (and thus uniformly), Ω ε are smooth for ε > 0 and the mean curvature of Ω ε is strictly positive.

We state below a few standard regularity results for regular domains. Proofs can be found, for example, in Theorems 2.2.2.3 and 3.1.2.1 from [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF].

Proposition 2.2. 1. If Ω is a domain with C 2 boundary. Then for every f ∈ L 2 (Ω) the problem -∆u = f in Ω u = 0 on ∂Ω has a unique solution solution u ∈ H 2 (Ω). 2.
If Ω is an open, bounded subset of R 2 with C 2 boundary and non-negative mean curvature, then there exists a constant C(Ω), which depends only on the diameter of Ω, such that

u H 2 (Ω) ≤ C(Ω) ∆u L 2 (Ω) , for every u ∈ H 2 (Ω) ∩ H 1 0 (Ω).
Now we are ready to prove the desired regularity result.

Theorem 2.3. Suppose Ω is a C 1,α set with non-negative weak mean curvature in L 2 (∂Ω).

Then for every f ∈ L 2 (Ω) the problem

-∆u = f in Ω u = 0 on ∂Ω has a unique solution u ∈ H 2 (Ω).
Proof: Proposition 2.1 implies the existence of a sequence of C 2 sets Ω n such that Ω n ⊂ Ω, Ω n has non-negative mean curvature and d H (∂Ω n , ∂Ω) → 0. We consider the solution

u m ∈ H 2 (Ω m ) of the Dirichlet problem in Ω m : -∆u m = f in Ω m u m = 0 on ∂Ω m ,
solution which exists by Proposition 2.2. Applying Theorem 1.5.1.5 from [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF] we deduce that ũm ∈ H 1 (R d ), where ũm is the extension by zero of u m outside Ω m . Proposition 2.2 implies that there exists a constant C such that u m H 2 (Ωm) ≤ C. This implies that ũm is bounded in

H 1 (R d ) and that the sequences v m,i,j = ∂ i ∂ j u m are bounded in L 2 (R d ).
Therefore, up to choosing a subsequence, we can assume that there exist

U ∈ H 1 (R d ) and V i,j ∈ L 2 (R d ) such that u m U weakly in H 1 (R d ) and v m,i,j V i,j weakly in L 2 (R d ).
In the following we denote by u the restriction of U to Ω. Since the supports of ũm are all contained in Ω, it follows that U is also supported on Ω. Thus u = 0 on ∂Ω. Now let ϕ ∈ C ∞ c (Ω) be a smooth function with compact support in Ω. Then there exists m 0 such that for every m ≥ m 0 the support of ϕ is contained in Ω m . Thus, for all m ≥ k 0 we have

Ω f ϕ = - Ωm ∆u m ϕ = Ωm ∇u m • ∇ϕ = Ω ∇ũ m • ∇ϕ.
Taking the limit as m → ∞ we get

Ω f ϕ = Ω ∇u • ∇ϕ.
Since this is true for every ϕ ∈ C ∞ c (Ω) we conclude that -∆u = f on Ω in the sense of distributions. Until now we have the existence of a solution u ∈ H 1 (Ω). In order to complete the proof we need to prove that the second derivatives of u are in L 2 . We take again ϕ ∈ C ∞ c (Ω) and for m ≥ m 0 we have

Ω ũm ∂ i ∂ j ϕ = Ωm u m ∂ i ∂ j ϕ = Ωm ∂ i ∂ j u m ϕ = Ω v m,i,j ϕ.
Taking the limit as m → ∞ we obtain

Ω u∂ i ∂ j ϕ = Ω V i,j ϕ,
which means that the distributional derivatives of order 2 denoted by ∂ i ∂ j u are given by V i,j and they are in L 2 for all i, j = 1...d. Thus we have proved the existence of a solution u ∈ H 2 (Ω) ∩ H 1 0 (Ω). The uniqueness follows at once using the estimates in Proposition 2.2 along with an approximation argument.

We note that the above arguments work if Ω does not have a singular part Σ. Indeed, the result stated in Proposition 2.1 is for the case when Ω is globally C 1 . When we have a singular set Σ it is no longer possible to use this. We can however prove that outside the singularity the eigenfunction is locally H 2 . In order to do this, let's recall the following integration by parts formula.

Ω ϕ|D 2 u| 2 dx + ∂Ω ϕH|∇u| 2 dσ = Ω ϕ(∆u) 2 dx + Ω (∇ϕ • ∇u)∆udx - Ω ∇ϕ • D 2 u∇udx.
For the sake of completness an idea of the proof, as well as some references are given in the Appendix C. Note that for a solution of (1) we know that the corresponding eigenfunction u is Lipschitz so ∇u ∈ L ∞ (Ω). See [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] for details. Moreover, since -∆u = λ k u and u ∈ L ∞ (Ω) we may as well assume that ∆u ∈ L ∞ (Ω). See for example [START_REF] Davies | Heat kernels and spectral theory[END_REF]Example 2.1.8] where we have the estimate u ∞ ≤ e 1/8π λ k (Ω) d/4 . Let x 0 be a point outside the singular set Σ. Since this set is closed, there is an open neighborhood of x 0 which does not intersect the singular part. Let B r be a ball centered in x 0 of radius r small enough such that B r ∩ ∂Ω is C 1,α . Consider a smooth cutoff function 0 ≤ ϕ ≤ 1 such that, ϕ = 1 on B r/2 and ϕ = 0 on ∂B r , where B r/2 is a ball concentric with B r , having radius r/2. With these considerations note that ∂(Ω∩Br) ϕH|∇u| 2 dσ ≥ 0 since on ∂Ω ∩ B r the curvature is non-negative in the distributional sense and on ∂B r ∩ Ω we have ϕ = 0. Thus, applying the above integration by parts formula for Ω ∩ B r and ϕ described above, we get

Ω∩Br ϕ|D 2 u| 2 dx ≤ Ω∩Br ϕ(∆u) 2 dx + Ω∩Br (∇ϕ • ∇u)∆udx - Ω∩Br ∇ϕ • D 2 u∇udx.
Applying the inequality ab ≤ a 2 +b 2 2 and Cauchy-Schwarz we obtain

Ω∩Br ϕ|D 2 u| 2 dx ≤ Ω∩Br ϕ(∆u) 2 dx + 1 2 Ω∩Br (∆u) 2 dx + 1 2 Ω∩Br (∇ϕ • ∇u) 2 dx + Ω∩Br ∇ϕ • D 2 u∇udx.
Integrating again by parts the last term in the above inequality, we obtain

Ω ∇ϕ • D 2 u∇udx = n i,j=1 Ω ∂ i ϕ∂ i ∂ j u∂ j u ∂ j = n i,j=1 ∂Ω ∂ i u∂ j u∂ i ϕν j dσ - n i,j=1 Ω ∂ i u∂ i ∂ j ϕ∂ j udx - n i,j=1 Ω ∂ i u∂ i ϕ∂ 2 j udx = ∂Ω (∇u • ∇ϕ)(∇u • n)dσ - Ω ∇u • D 2 ϕ∇udx - Ω (∇u • ∇ϕ)∆udx.
Since ∇u is in L ∞ (Ω) and in our case -∆u = λ k u ∈ L ∞ (Ω), the above expression is bounded and we may conclude that Ω∩Br ϕ|D 2 u| 2 dx < +∞ which means that u ∈ H 2 (B r/2 ). Since this argument is valid for any point x 0 outside Σ, it follows that u ∈ H 2 loc (Ω \ Σ).

Analytic perturbations for eigenvalue problems

We know that the eigenvalues of the Dirichlet Laplace operator on Ω are differentiable if they are simple. We refer to [START_REF] Henrot | Variation et optimisation de formes[END_REF] for further details. However, numerical results found in [START_REF] Antunes | Optimisation of eigenvalues of the Dirichlet Laplacian with a surface area restriction[END_REF][START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF] show that optimizers of problem (1) tend to have multiple eigenvalues. It is possible to recover some information regarding the differentiability of the eigenvalues, even when the multiplicity is greater than one, by using the theory of operator perturbations.

Rellich proved in [START_REF] Rellich | Störungstheorie der Spektralzerlegung[END_REF] that analytic perturbations of a self-adjoint operator with the same domain allows an analytic parametrization of the eigenvalues and the eigenfunctions. More recent results can be found in [START_REF] Kriegl | Differentiable perturbation of unbounded operators[END_REF]. The article of Micheletti [START_REF] Micheletti | Perturbazione dello spettro dell'operatore di Laplace, in relazione ad una variazione del campo[END_REF] contains details about the application of Rellich's result to the study of eigenvalue problems. Indeed, it is proved that the Laplace operator defined on Ω ε = (Id + εV )(Ω) has the same spectrum as a self-adjoint operator A ε which depends analytically on ε and has a domain of definition independent of ε. The result of Rellich, stated above, implies that the eigenvalues and eigenvectors of (-∆) on Ω ε can be parametrized analytically with respect to ε. Thus, if λ is an eigenvalue of Ω with multiplicity m, with the associated orthonormal basis (u i ) m i=1 , there exist m functions ε → λ ε,i and the associated family of eigenfunctions ε → ϕ ε,i orthonormal in L 2 (Ω ε ), such that both dependencies are analytic in ε, λ 0,i = λ, ϕ 0,i = u i and

-∆ϕ ε,i = λ ε,i ϕ ε,i in Ω ε ϕ ε,i = 0 on ∂Ω ε ,
for every i = 1, ..., m. Differentiating with respect to ε we obtain the classical derivative formula

dλ ε,i dε ε=0 = - ∂Ω ∂u i ∂n 2 V.ndσ (3)
For details see [START_REF] Oudet | Numerical minimization of eigenmodes of a membrane with respect to the domain[END_REF] where a formal argument shows how to deduce relation [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF]. The formal argument is rigorous once we know that the eigenvalues and eigenfunctions depend analytically on ε. The H 2 regularity (local H 2 when singularities are present) of the eigenfunctions proved in the previous section allows writing the derivative of the eigenvalues using the above boundary integral depending on the normal derivative. For details see [START_REF] Henrot | Variation et optimisation de formes[END_REF]Theorem 5.7.4]. Proof: Since Ω * has C 1,α regularity outside an eventual closed singular set Σ of dimension at most d -8, it is possible to represent the boundary around a point x 0 ∈ ∂Ω * \ Σ as the graph of a C 1,α function h : B a ⊂ R d-1 → R, where B a is a ball of radius a, small enough such that we are still outside the singular set. We denote by Ω ε = (Id + εV )(Ω * ) for a C 1 (R d , R d ) vector field. Using this type of graph representation we deduce that if ψ is the perturbation of this graph associated to V then ε → Per(Ω ε ) is differentiable at zero and its derivative is equal to

Main result

Ba ∇ψ • ∇h 1 + |∇h| 2 .
As noted in the previous section, the function ε → λ k (Ω ε ) has left and right derivatives. Furthermore, we know that there is an orthonormal basis of the associate eigenspace, denoted (u i ) p i=1 , such that the directional derivatives of λ k are in the set

- Ba |∇u i | 2 ψ(x)dx p i=1
. This is just a simple consequence of the derivative formula (3) in the local coordinate system. Indeed, the normal derivative is just the gradient since the tangential derivative is zero. On the other hand the term V.n is equal to ψ(x)/ 1 + |∇h| 2 and the Jacobian is 1 + |∇h| 2 .

Thus, at the optimum, regardless of the multiplicity and the C 1 perturbation considered, there exist left and right derivatives for λ k (Ω ε )+Per(Ω ε ) at ε = 0. The case where the eigenvalue is simple is straightforward and is similar to the approach used in [START_REF] Bucur | Minimization of λ 2 (Ω) with a perimeter constraint[END_REF]. Suppose now that λ k is multiple. The local optimality of Ω * implies that the left and right derivatives of λ k (Ω ε )+Per(Ω ε ) at ε = 0 in any direction given by ψ in the local coordinates are of different signs. As a consequence, given a perturbation ψ ∈ C 1 (B a ) we have two indices i, j such that Ba ∇ψ • ∇h

1 + |∇h| 2 - Ba |∇u i | 2 ψ ≥ 0 (4) Ba ∇ψ • ∇h 1 + |∇h| 2 - Ba |∇u j | 2 ψ ≤ 0 (5)
We define the linear functionals F, G i on C 1 (B a ) by

F(ψ) = Ba ∇ψ • ∇h 1 + |∇h| 2 G i (ψ) = Ba |∇u i | 2 ψ
Equations ( 4) and ( 5) tell us that for every ψ ∈ C 1 (B a ) there exist indices i, j such that F(ψ) ≥ G i (ψ) and F(ψ) ≤ G j (ψ). This implies, in particular, that for every ψ ∈ C

1 + |∇h| 2 = m i=1 µ i Ba |∇u i | 2 ψ,
for every ψ ∈ C 1 (B a ). Since h is already C 1,α , we know that each u i is also C 1,α [13, Theorem 6.18]. Therefore ∇u i are all in C 0,α .

In [START_REF] Landais | A regularity result in a shape optimization problem with perimeter[END_REF] it is proved that if h ∈ C 1,α (B a ). f ∈ C 0,α (B a ) and the distributional equation

Ba ∇ψ • ∇h 1 + |∇h| 2 = Ba f ψ, ∀ψ ∈ C 1 (B a )
holds, then h ∈ C 2,α (B a ). Using this result we conclude that h ∈ C 2,α and h is a strong solution of the equation

-div ∇h 1 + |∇h| 2 = m i=1 µ i |∇u i | 2 (6)
with µ 1 + ... + µ p = 1. Using straightforward computations equation ( 6) is equivalent to

- ∆h 1 + |∇h| 2 + ∇h • ∇h (1 + |∇h| 2 ) 3/2 = m i=1 µ i |∇u i | 2
which leads to

-∆h = 1 + |∇h| 2 m i=1 µ i |∇u i | 2 - |∇h| 2 1 + |∇h| 2 . ( 7
)
Since h ∈ C 2,α (B a ) and this is true for any local chart which does not intersect Σ, the set Ω is a domain with C 2,α boundary outside Σ. Therefore using again the standard Schauder regularity results [START_REF] Gilbarg | Elliptic partial differential equations of second order[END_REF]Theorem 6.18] we see that the eigenfunctions (u i ) p i=1 are in C 2,α (B a ). This means that the right hand side of equation ( 7) is in C 1,α (B a ). Theorem 9.19 from [START_REF] Gilbarg | Elliptic partial differential equations of second order[END_REF] allows us to deduce that h belongs to C 3,α (B a ). In general we see that if h ∈ C k,α (B a ), k ≥ 2, then the right hand side of ( 7) is in C k-1,α (B a ) and thus h ∈ C k+1,α (B a ). An inductive bootstrap argument allows us to see that h ∈ C ∞ (B a ). Moreover, since the coefficients of the partial differential equation ( 7) are analytic and h is C ∞ the results of Morrey [START_REF] Morrey | On the analyticity of the solutions of analytic non-linear elliptic systems of partial differential equations. I. Analyticity in the interior[END_REF] allow us to deduce that h is analytic. Repeating the argument around each point outside Σ, we conclude that ∂Ω \ Σ is analytic.

Remark 4.2. Theorem 4.1 improves the results in [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF] in yet another direction. In the cited reference the authors prove that at the optimum there exists a family of eigenfunctions (ϕ i ) m i=1 in the eigenspace associated to λ k such that

H = m i=1 (∂ n ϕ i ) 2 .
In the proof of our result we deduce a more precise result which says more about the family (ϕ i ) and about the number of eigenfunctions present in the optimality condition. Indeed we prove that if (u i ) m i=1 is an orthonormal basis of the eigenspace associated to λ k then there exist µ 1 , ..., µ m ∈ [0, 1] with µ 1 + ... + µ m = 1 such that

H = m i=1 µ i (∂ n u i ) 2 .
We note that this behavior has been conjectured by the numerical results presented in [START_REF] Bourdin | Optimal partitions for eigenvalues[END_REF]. (P2) F is of class C 1 and at least one of its partial derivatives does not vanish when evaluated at (λ k 1 (Ω * ), ..., λ kp (Ω * )).

(P3) F is increasing in each variable, furthermore, for any compact K ⊂ R p \ {0} there exists a > 0 such that if x, y ∈ R p with x j ≥ y j , j = 1, ..., p then F (x) -F (y) ≥ a|x -y|.

Properties (P1) and (P3) are taken from [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] to guarantee existence. Property (P2) is stronger than the analogous property in [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] (locally Lipschitz continuity) and allows us to differentiate the functional at the optimum.

Remark 4.4. The techniques presented here should also apply to the problem min{λ k (Ω) + Per(Ω), |Ω| = m} introduced in [START_REF] De Philippis | Regularity of minimizers of shape optimization problems involving perimeter[END_REF]. Indeed, the functional is the same, but there is a volume constraint. The existence of a solution which is C 1,α regular outside a closed singular set of dimension at most d -8 has been treated in [START_REF] De Philippis | Regularity of minimizers of shape optimization problems involving perimeter[END_REF]. Like above it is possible to write the expression of the shape derivative, but we need to work with flows which preserve the volume. Such flows can be constructed for any perturbation V such that ∂Ω V.n = 0 (see [START_REF] Soufi | Domain deformations and eigenvalues of the Dirichlet Laplacian in a Riemannian manifold[END_REF]). Then, using techniques similar to the ones in [START_REF] Soufi | Domain deformations and eigenvalues of the Dirichlet Laplacian in a Riemannian manifold[END_REF] we may deduce that at the optimality condition (6) becomes

-div ∇h 1 + |∇h| 2 = m i=1 µ i |∇u i | 2 + const.
This allows us to deduce the full regularity of the above problem in the same way as in Theorem 4.1.

The results in this section allow us to give a different argument for [6, Theorem 2.5] and extend the corresponding result to all dimensions. The following result also generalises [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF]Theorem 5.7]. Proof: This is a simple consequence of the fact that ∂Ω * is analytic outside a singular set of dimension at most d -8.

The Multiplicity Cluster

Once the regularity of the minimizers of ( 1) is proved we may look into other facts about the optimal set. As noted in [START_REF] Antunes | Optimisation of eigenvalues of the Dirichlet Laplacian with a surface area restriction[END_REF] and [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF] the numerical experiments show that in most cases the multiplicity of the k-th eigenvalue of the optimal set is larger than 1. Numerical experiments also suggest that when λ k (Ω * ) is multiple then λ k (Ω * ) < λ k+1 (Ω * ). A first result in this sense has been obtained in [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF] and is similar to [15, Lemma 2.5.9]. For the sake of completness we give a sketch of the proof.

Theorem 5.1. Suppose Ω * is a minimizer for problem (1). If we have λ k-1 (Ω * ) < λ k (Ω * ) then the k-th eigenvalue is simple, i.e. λ k (Ω * ) < λ k+1 (Ω * ).

Proof: Suppose λ k (Ω * ) is multiple and λ k (Ω * ) = ... = λ k+m-1 (Ω * ). Since λ k is the smallest eigenvalue in the cluster and Ω * is a local minimizer, the m analytic parametrizations of the eigenvalues of the spectrum, described in Section 3, have derivative zero at ε = 0. This would imply that (∂ n u i ) 2 = H for all i = 1, ..., m.

If m > 1 this implies the existence of two different eigenfunctions u, v associated to

λ k (Ω * ) such that (∂ n u) 2 = (∂ n v) 2 = H.
The Hölmgren uniqueness theorem implies that H cannot vanish on a relatively open set of ∂Ω (see [START_REF] Bogosel | Qualitative and numerical analysis of a spectral problem with perimeter constraint[END_REF]Theorem 5.7] for an alternative argument). Thus, there exists a relatively open set γ ⊂ ∂Ω such that H > 0 on γ. We can see that in this situation

∂ n u = ± √ H and ∂ n v = ± √ H on γ.
This implies that we have either ∂ n (u + v) = 0 or ∂ n (u -v) = 0 on γ. By Hölmgren's uniqueness theorem we deduce that u = ±v which is a contradiction. Therefore λ k (Ω * ) is simple.

We can provide a stronger result based on the methods used in Theorem 2.5.8 and Lemma 2.5.9 from [START_REF] Henrot | Extremum problems for eigenvalues of elliptic operators[END_REF]. We note that we have already proved that optimal sets Ω * are analytic outside a singular set Σ so the cited results apply in our case. Below we prove a generalization of Lemma 2.5.9 from [START_REF] Henrot | Extremum problems for eigenvalues of elliptic operators[END_REF]. Lemma 5.2. Suppose Ω is a smooth set in R d and that λ k+1 (Ω) = ... = λ k+m (Ω). Then there exist vector fields V such that of the derivatives dλ k+i (Ω ε ) dε ε=0 are strictly negative and the other m -are strictly positive, for ∈ {1, 2, ..., m -1}. Moreover, V can be chosen in such a way that the derivative of the perimeter of Ω in the direction of V is zero, which means that ∂Ω HV.ndσ = 0.

Remark 5.3. Working like in [3, Lemma 2.6] we may construct a family of sets Ω ε = φ ε (Ω) such that Ω ε have the same perimeter as Ω = Ω 0 and the derivative of the flow φ ε at ε = 0 is V .

Proof of Lemma 5.2: Theorem 2.5.8 from [START_REF] Henrot | Extremum problems for eigenvalues of elliptic operators[END_REF] says that the derivatives of the eigenvalues of the perturbations Ω ε of Ω are the eigenvalues of the m × m matrix

A V = (a i,j ) where a i,j = - ∂Ω ∂ n u i ∂ n u j V.ndσ, 1 ≤ i, j ≤ m.
As a consequence, the derivatives of λ k (Ω * ) + Per(Ω * ) in the direction V are the eigenvalues of A shifted by the derivative of the perimeter for the same direction V .

As in the proof of the previous theorem we place ourselves in a relatively open subset γ of ∂Ω where the mean curvature H is strictly positive. For a point X ∈ γ we denote ψ(X) = (∂ n u 1 (X), ..., ∂ n u m (X)) ∈ R m . We can find points X 1 , ..., X n ∈ γ such that the vectors ψ(X i ), i = 1, ..., m are linearly independent. We present a simple proof of this fact in Appendix A.

Once chosen X 1 , ..., X n we consider a vector field V δ such that V δ is supported on O 1 ∪...∪O n , where X p ∈ O p and H n-1 (O p ) = δ. Given j ∈ {1, 2, ..., n} we can consider V δ to be of the form

V δ = (k p /δ) n(x) on O p \ Q p (f p (x)/δ) n(x) on Q p
where constants k p satisfy k p < 0 for p ≤ and k p > 0 for p > , Q p are portions around the boundary of O p of H d-1 measure O(δ 2 ), f p smooth extensions of the constant such that f p ≤ k p and n(x) denotes the normal to ∂Ω. Furthermore, the constants k i are chosen such that ∂Ω HV δ .ndσ = 0, so that the derivative of the perimeter is zero. This can always be done. It suffices to fix the first constants k i = c < 0 and choose the rest equal to a positive constant which makes the above integral zero. Note that the derivative of the perimeter cannot be equal to zero if we choose all k i to have the same sign. This is due to the fact that γ is a region where H is strictly positive.

It is not difficult to see that

lim δ→0 Op ∂ n u i ∂ n u j V δ .ndσ = k p ∂ n u i (X i )∂ n u j (X i ), p = 1, ..., n.
Therefore, the matrix A V δ converges to

k 1 ψ(X i ) T ψ(X i ) + ... + k n ψ(X n ) T ψ(X n )
as δ goes to 0. For a general vector X ∈ R n we have

X T A V δ X → k 1 (X.ψ(X 1 )) 2 + ... + k n (X.ψ(X n )) 2 .
Since the vectors ψ(X i ) are linearly independent we deduce that this limit matrix is associated to a quadratic form of signature (n -, ) so it has negative eigenvalues and n -positive ones. The continuity of the eigenvalues of the matrices with respect to the entries implies that for δ small enough A V δ has negative eigenvalues and n -positive ones. Thus eigenvalues from the cluster have negative derivatives and n -have positive derivatives. Note that even if the above result is stated for smooth sets, we may apply it in our case by performing only perturbations supported outside the singularity. An immediate corollary is the following.

Corollary 5.4. If Ω * is an optimizer for (1) corresponding to the index k then λ k (Ω * ) < λ k+1 (Ω * ).

Proof: Suppose that λ k (Ω * ) is multiple and λ s (Ω * ) = ... = λ k (Ω * ) = λ k+1 (Ω * ) = ... = λ S (Ω * ). Apply Lemma 5.2 to deduce that it is possible to find a perturbation V of Ω * such that the derivatives of the eigenvalues in the cluster are strictly negative for indices in [s, k] and strictly positive for indices in [k + 1, S]. Moreover, the perturbation can be chosen so that the perimeter has derivative zero in the direction of V . This perturbation would then decrease λ k + Per, thus contradicting optimality. Therefore we can only have λ k (Ω * ) < λ k+1 (Ω * ). Denote Ω k and Ω k+1 solutions of problem [START_REF] Ambrosio | Functions of bounded variation and free discontinuity problems[END_REF] for k and k + 1, respectively. Then we have the following inequalities

λ k (Ω k ) + Per(Ω k ) ≤ λ k (Ω k+1 ) + Per(Ω k+1 ) ≤ λ k+1 (Ω k+1 ) + Per(Ω k+1 ),
where the first inequality comes from the optimality of Ω k and the second from the ordering of the eigenvalues. If the optimal costs c k , c k+1 defined in (2) satisfy c k = c k+1 then the above inequalities imply that λ k (Ω k+1 ) = λ k+1 (Ω k+1 ) and that Ω k+1 is also solution of (1) for k. This is in contradiction with the previous corollary.

Remark 5.6. Remark 5.3 allows us to conclude a similar result for the constrained problem

d k = min{λ k (Ω) : Per(Ω) = c}.
In this case we also have d k < d k+1 .

Remark 5.7. Using the results above it is possible to deduce that for the shape Ω * which minimizes λ 2 (Ω) + Per(Ω) the second eigenvalue is simple and thus we have (∂ n u 2 ) 2 = H. Indeed, the optimal set Ω * is connected [START_REF] Philippis | Existence and regularity of minimizers for some spectral functionals with perimeter constraint[END_REF] so we have λ 1 (Ω * ) < λ 2 (Ω * ). The results above show that λ 2 (Ω * ) < λ 3 (Ω * ) so λ 2 (Ω * ) is simple.

A From discrete to continous linear dependence Consider a positive integer n ≥ 2 and let f 1 , f 2 , ..., f n be real continuous functions defined on an open connected set I ⊂ R d . For a point x ∈ I we denote ψ(x) = (f 1 (x), ..., f n (x)) ∈ R n . We suppose that none of the functions f i is identically zero on I. Suppose that for each set of n different points x 1 , ..., x n ∈ I the vectors ψ(x 1 ), ..., ψ(x n ) are linearly dependent on an open connected subset of I. Then the functions f 1 , ..., f n are linearly dependent, i.e. there exist constants α 1 , ..., α n not all zero such that

α 1 f 1 + ... + α n f n = 0.
Let's start with the case n = 2. Pick x 1 = x 2 in I. We know that the vectors ψ(x 1 ), ψ(x 2 ) are linearly dependent which means that

f 1 (x 1 ) f 2 (x 1 ) = λ f 1 (x 2 ) f 2 (x 2 )
We find that on the set {f 2 = 0} the function f 1 /f 2 is constant. Thus, picking eventually a connected component of {f 2 = 0} we find that f 1 , f 2 are linearly dependent on a connected open set.

We can prove the result for general n ≥ 2 by induction. Indeed, let's suppose that the result holds for n functions. Consider now f 1 , ..., f n+1 defined on I such that for every n + 1 points x 1 , ..., x n+1 the vectors ψ(x 1 ), ..., ψ(x n+1 ) are linearly dependent. Suppose that no two functions f i are equal or else we have nothing to prove. Furthermore, we choose x 1 such that f (x 1 ) = 0. Thus there exist scalars λ 1 , ..., λ n , not all zero such that

           λ 1 f 1 (x 1 ) + λ 2 f 1 (x 2 ) + ... + λ n+1 f 1 (x n+1 ) = 0 λ 1 f 2 (x 1 ) + λ 2 f 2 (x 2 ) + ... + λ n+1 f 2 (x n+1 ) = 0 • • • λ 1 f n+1 (x 1 ) + λ 2 f n+1 (x 2 ) + ... + λ n+1 f n+1 (x n+1 ) = 0.
Without loss of generality suppose that λ 1 = 0. Eliminating the elements below λ 1 f 1 (x 1 ) in the first column leaves us with a system of the form

     µ 2 (f 2 -f 1 )(x 2 ) + ... + µ n+1 (f 2 -f 1 )(x n+1 ) = 0 • • • = 0 µ 2 (f n+1 -f 1 )(x 2 ) + ... + µ n+1 (f n+1 -f 1 )(x n+1 ) = 0
Therefore f 2 -f 1 , ..., f n+1 -f 1 satisfy the induction hypothesis and thus they are linearly dependent. It is obvious that this implies that f 1 , ..., f n+1 are linearly dependent.

B A variant of the Hahn-Banach Theorem

We recall below the well known Hahn-Banach separation theorem. For a proof or more details see [START_REF] Clarke | Functional analysis, calculus of variations and optimal control[END_REF]Section 2.4].

Theorem B.1. Ket K 1 , K 2 be nonempty, disjoint convex subsets of the normed vector space X.

Then if K 1 is open, there exist ζ ∈ X * and θ ∈ R such that ζ, x < θ ≤ ζ, y ∀x ∈ K 1 , y ∈ K 2 .
An immediate consequence is the following. Proposition B.2. Let {ζ i : i = 1, 2, ..., k} be a finite subset in X * . The following are equivalent: (a) There is no ν ∈ X such that ζ i , ν < 0, ∀i ∈ {1, 2, ..., k}.

(b) The set {ζ i : i = 1, 2, ..., k} is positively linearly dependent: there exists a non zero nonnegative vector γ ∈ R k such that k i=1 γ i ζ i = 0.

Proof: First let's note that we may suppose that none of the functionals ζ i are identically zero, since the problem would reduce itself to a smaller value of k or to the trivial case when all functionals considered vanish.

We start with the implication (b) =⇒ (a). Suppose that k i=1 γ

i ζ i = 0. Now, if ν ∈ X we have 0 = k i=1 γ i ζ i , ν ≤ ( k i=1 γ i ) max i ζ i , ν , which implies max i ζ i , ν ≥ 0, since k i=1 γ i > 0.
For the implication (a) =⇒ (b) we consider the following subsets of R k :

K 1 = {y ∈ R k : y i < 0∀i ∈ {1, 2, ..., k}} K 2 = {( ζ 1 , ν , ζ 2 , ν , ..., ζ k , ν ) : ν ∈ X}.
We see immediately that (a) implies K 1 ∩ K 2 = ∅ and furthermore K 1 is open. Therefore, we can apply the Hahn-Banach separation theorem B.1 and find a functional ϕ : R k → R and θ ∈ R such that ϕ, x < θ ≤ ϕ, y , ∀x ∈ K 1 , y ∈ K 2 .

We know that such a functional is of the form ϕ(x) = Taking y n i = -n, y n j = -1/n, j = i and taking n → ∞ we see that this is possible only if all γ i are non-negative. This finishes the proof.

C A formula relating the Laplacian and the Hessian Theorem 3.1.1.1. from [START_REF] Grisvard | Elliptic problems in nonsmooth domains[END_REF] says that for v ∈ H 1 (Ω) n and Ω of class C 2 we have

Ω | div v| 2 dx - n i,j=1 Ω ∂ i v j ∂ j v i dx = (tangential part) + ∂Ω H(v • nν) 2 dσ,
where H is the mean curvature of Ω. If we put v = ∇u for u = 0 on ∂Ω, then the tangential part is zero and we obtain a relation between the Laplacian and the Hessian. The aim here is to write a similar relation when we have an additional multiplication with a function ϕ ∈ H 1 (R d ).

We wish to prove the following The terms integrated on Ω are

A Ω = n i,j=1 Ω ∂ j ϕv i ∂ i v j dx - n i,j=1 Ω v i ∂ i ϕ∂ j v j dx = Ω {(v • ∇)v • ∇ϕ}dx - Ω (v • ∇ϕ) div vdx.
In order to avoid possible complications, unnecessary in our case, we recall that we want to use this computation for v = ∇u with u = 0 on ∂Ω. This means, in particular, that the tangential part of v is v τ = (∇u) τ = 0. Therefore, in the following, we suppose that the tangential part of v is zero. Using the reasoning from [14, Theorem 3.1.1.1] and neglecting the tangential components we obtain Note that even if the result above is stated for C 2 sets, we may apply it in the case where Ω is C 1,α with non-negative and bounded distributional curvature. Analysing the proof we notice that the only difficulty comes from the boundary integrals A ∂Ω where which are considered for fields v with zero tangential components.

A ∂Ω = -

  min |Ω|=m (Per(Ω) + G(Ω))

Theorem 4 . 1

 41 If Ω * is a local minimizer for problem (1) then Ω * is analytic for d ≤ 7. If d ≥ 8 then Ω * is analytic outside a closed set of Hausdorff dimension d -8.

Theorem 4 . 1 .

 41 If Ω * is a local minimizer for problem (1) then Ω * is analytic for d ≤ 7. If d ≥ 8 then Ω * is analytic outside a closed set of Hausdorff dimension at most d -8.

Remark 4 . 3 .

 43 The results of Theorem 4.1 can be generalized to problems of the type minΩ∈R d F (λ k 1 (Ω), ..., λ kp (Ω)) + Per(Ω)where F : R p → R + satisfies the properties (P1) F (x) → ∞ as |x| → ∞.

Corollary 4 . 5 .

 45 If Ω * is a solution for problem (1) then • ∂Ω * does not contain flat parts • ∂Ω * does not contain open subsets of spheres unless Ω * is a ball.

Remark 5 . 5 .

 55 A direct consequence of the previous theorem is that the optimal cost for problem (1) is a strictly increasing function of k, i.e. min Ω∈R d λ k (Ω) + Per(Ω) < min Ω∈R d λ k+1 (Ω) + Per(Ω).

γ

  i x i , where γ = (γ 1 , ..., γ k ) is a non-zero element of R k . Therefore ϕ, y ≥ θ, ∀y ∈ K 2 becomes θ ≤ k i=1 γ i ζ i , ν , ∀ν ∈ X,which is only possible when k i=1 γ i ζ i = 0. On the other hand, we have k i=1 γ i y i < θ for all y = (y i ) ∈ K 1 .

Proposition C. 1 .

 1 If Ω is of class C 2 , ϕ ∈ H 1 (R d ) and u is smooth in Ω ∩ supp ϕ then Ω ϕ|D 2 u| 2 dx + ∂Ω ϕH|∇u| 2 dσ = Ω ϕ(∆u) 2 dx + Ω (∇u • ∇ϕ)∆udx -Ω ∇ϕ • D 2 uDudxProof: We integrate by parts two times and we obtainΩ ϕ| div v| 2 = n i,j=1 Ω ϕ∂ i v i ∂ j v j dx Ω v i ϕ∂ i,j v j dx -n i,j=1 Ω v i ∂ i ϕ∂ j v j dx + n i,j=1 ∂Ω ϕv i ∂ j v j ν i dσ ∂ j = n i,j=1 Ω ϕ∂ j v i ∂ i v j dx + n i,j=1 Ω ∂ j ϕv i ∂ i v j dx -n i,j=1 ∂Ω ϕv i ∂ i v j ν j dσ--n i,j=1 Ω v i ∂ i ϕ∂ j v j dx + n i,j=1 ∂Ω ϕv i ∂ j v j ν i dσ Therefore we have I(v) = Ω ϕ| div v| 2 -n i,j=1 Ω ϕ∂ j v i ∂ i v j dx= A Ω (terms on Ω) + A ∂Ω (terms on ∂Ω)

∂Ω ϕHv 2 dσ

 2 ϕv i ∂ i v j ν j dσ + n i,j=1 ∂Ω ϕv i ∂ j v j ν i dσ = ∂Ω ϕv ν div vdσ -∂Ω ϕ{(v • ∇)v} • νdσ =As stated before, we replace v by ∇u where u = 0 on ∂Ω and we getΩ ϕ(∆u) 2 dx -Ω ϕ|D 2 u| 2 dx = Ω ∇ϕ • D 2 uDudx -Ω (∇u • ∇ϕ)∆udx + ∂Ω ϕH|∇u| 2 dσ. Therefore Ω ϕ|D 2 u| 2 dx + ∂Ω ϕH|∇u| 2 dσ = Ω ϕ(∆u) 2 dx + Ω (∇u • ∇ϕ)∆udx -Ω ∇ϕ • D 2 uDudx

Acknowledgements. The author wishes to thank Dorin Bucur, Otis Codosh, Jimmy Lamboley and Bozhidar Velichkov for useful discussions in connection with this work. The author also thanks the anonymous reviewers for their careful analysis of the text and for the suggestions which helped improve the quality of this paper. This work was supported by the project ANR Optiform and by the Fondation Sciences Mathematiques de Paris.