Yannis Merlet

Roland Baviere

Nicolas Vasset

Validation of a multi-objective optimization framework for the sizing of pipes in DH Networks

Keywords: district heating network, multi-objective optimization, sizing, simulation District Heating GA: Genetic Algorithm CAPEX: Capital Expenditure OPEX: Operational Expenditure

Sizing distribution equipment of a district heating (DH) network is a complex yet decisive process to target techno-economic optimality for these systems. This paper describes and validates a framework that uses multi-objective optimization to support decisions regarding the sizing of a DH network. A genetic algorithm (NSGA2) is employed to generate optimal pipe diameters with respect to the operational and investment cost. This framework is then validated using an elementary toy problem consisting of a straight horizontal two-tube network, with a set of identical consumers equally distributed in the network and supplied by a unique heat source. We show that the framework, on the aforementioned problem, reaches very good convergence and diversity for the optimal solutions. This sizing method is then compared to typical local sizing method. We demonstrate a significant decrease of global cost by 40% by using DH network system simulation to take into account local interactions inside the DH system.

Introduction

District heating (DH) networks are an efficient way of producing and distribute heat on a territory and are widely used in Europe [START_REF] Werner | International review of district heating and cooling[END_REF]. They have the advantages of increased heat production efficiency and the possibility to distribute heat from various sources such as boilers, renewable energy, or heat recovery in the industry. However, the design of such DH network is a challenge on the production side, and on the distribution side. For example, oversizing the distribution pipes can lead to useless investment cost and under sizing may prevent some customer's demand to be fulfilled. In this domain, engineering guidelines generally rely on a pipe by pipe approach using threshold on either the fluid velocity, the regular pressure loss or both [START_REF] Frederiksen | District Heating and Cooling[END_REF]. Recommendations may vary depending on the local standards but the calculations are often performed on a static sizing point.

Researchers proposed different methods to use optimization in order to generate optimal designs for a distribution network. Most of the proposed methods aim to size the layout of the distribution network while the possible trajectories of potential pipes are already predetermined [START_REF] Egberts | Challenges in heat network design optimization[END_REF]. Multiple studies, referenced in [START_REF] Sameti | Optimization approaches in district heating and cooling thermal network[END_REF], address this problem by simplifying the network, in order to solve the sizing problem as a MILP problem. This approach leads to affordable computational effort only when complex DH system featuring multiple producers or storage units are discarded.

In order to deal with non-linear problems, heuristic methods have been used in recent research. [START_REF] Vesterlund | Design Optimization of a District Heating Network Expansion, a Case Study for the Town of Kiruna[END_REF]) used a Genetic Algorithm (GA) to size the piping and a boiler for an extension of Kiruna's DH network while lowering the CAPEX and OPEX of the network. [START_REF] Li | District Heating Network Design and Configuration Optimization with Genetic Algorithm[END_REF]) used a GA to minimize the net present cost of the distribution of a DH network with 10 consumers by optimizing both the topology and the sizing of the network. Those methods are efficient to deal with non-linearities and non-convexities in the evaluation function, and can deal with discrete decision variable. Thus, it makes them suitable for our problem.

Heuristic methods have been used for the sizing of large water distribution networks with 317 decision variables (pipe diameters) by [START_REF] Monsef | Comparison of evolutionary multi objective optimization algorithms in optimum design of water distribution network[END_REF], and proved to have a good scalability to tackle problems such as the sizing of a DH network at a city-scale. Nevertheless, the computational cost of these methods can be high as [START_REF] Monsef | Comparison of evolutionary multi objective optimization algorithms in optimum design of water distribution network[END_REF] had to call 500 000 evaluations of its system for 317 decision variables and [START_REF] Vesterlund | Design Optimization of a District Heating Network Expansion, a Case Study for the Town of Kiruna[END_REF] called 25 000 evaluation for 111 decision variables. There is therefore an interest in the use of multi-objective optimization if the evaluation function is not expensive to calculate. This paper reports on the first steps towards the design and validation of an optimization framework dedicated to the sizing of distribution pipes in a DH network. The main contribution lies in the combination of an efficient DH simulation tool with a multi-objective optimization method and its validation against a case study. The combination of these elements enable to identify the best design choices in a discrete decision space. Used in combination with a scalable thermal-hydraulic simulation engine, this approach allows dealing with the sizing of large-scale DH networks while using meta-heuristic methods. In the present work, we rely on DistrictLab-H, for the simulation of DH network models. More details about this simulation engine are given in the next section and in a paper in preparation [START_REF] Bavière | A step towards predictive simulation of large-scale heating grids: Implementation and validation of the DistrictLab-H platform[END_REF]. This paper is organized as follows. The next section describes the optimization framework with its structure and components. Section 3 introduces a toy problem used for validation purpose. In Section 4, the results of the optimization are presented and compared to a reference solution while section 5 compares the results obtained using a local sizing method widely used. The last section concludes our work.

Framework

The main design goal of our optimization framework is to combine a multi-objective optimization method, with a non-linear evaluation function brought by a DH network simulation engine. This section first describes the specification of the problem; secondly, it presents the choices made in the framework and then explains the global structure of our framework.

Specifications of the optimization problem

The decision variables to optimize are tunable parameters belonging to a DH Network simulation model. In the present work, we only considered the distribution diameters but more parameters will be included in further research, such as the possibility to add new boilers, storage units or booster pumps. The targeted final application for the framework is the optimization of the DH network model for a medium city comprising several hundreds of customers and thousands of distribution pipes. Those parameters are by nature discrete, because not all diameters of distribution pipes exists in the market.

The objective functions presently optimized in the framework are the investment cost (CAPEX) for the proposed modification of the network and the operation cost (OPEX) of the network. However, additional objectives, like CO2 emissions for example, could also be considered.

For the optimization results to be applicable, some constraints have to be accounted for: first of all, a constraint for the satisfaction of all customers is integrated : this constraint is expressed under the form of a user defined satisfaction rate threshold:

𝑆𝑎𝑡𝑖𝑠𝑓𝑎𝑐𝑡𝑖𝑜𝑛 𝑟𝑎𝑡𝑒 = 1 𝑁 𝑐𝑙𝑖𝑒𝑛𝑡𝑠 ∑ 𝐻𝑒𝑎𝑡 𝑑𝑒𝑙𝑖𝑣𝑒𝑟𝑒𝑑 𝐻𝑒𝑎𝑡 𝑑𝑒𝑚𝑎𝑛𝑑 𝑐𝑙𝑖𝑒𝑛𝑡𝑠 ≥ 𝑆𝑎𝑡𝑖𝑠𝑓𝑎𝑐𝑡𝑖𝑜𝑛 𝑟𝑎𝑡𝑒 𝑡ℎ𝑟𝑒𝑠ℎ𝑜𝑙𝑑 1 2.2.

Evaluation function

The evaluation of a possible design solution is done through the DistrictLab-H DH network simulation tool. The corresponding model and resolution method are described in a paper in preparation [START_REF] Bavière | A step towards predictive simulation of large-scale heating grids: Implementation and validation of the DistrictLab-H platform[END_REF]. This tool allows the dynamic simulation of distribution pipes connected to centralized supply units and to a collection of heat consumers through substations of the indirect type. This solver enables the simulation of large-scale DH networks in a reasonable calculation cost.

Optimization method

The discrete nature of the decision parameters, the possible number of these decision parameters and the nonlinear nature of the evaluation function drive us towards heuristic multi-objective optimization methods. In those methods, GA are widely used for problems with the same characteristics in the literature [START_REF] Sameti | Optimization approaches in district heating and cooling thermal network[END_REF][START_REF] Vesterlund | Design Optimization of a District Heating Network Expansion, a Case Study for the Town of Kiruna[END_REF].

We chose to used NSGA-II genetic algorithm [START_REF] Deb | A fast and elitist multiobjective genetic algorithm: NSGA-II[END_REF] for its robustness and its performance on solving multi-objective optimization problems with a limited amount of objectives. Among the implementation of NSGA-II, the one from python library Deap [START_REF] Fortin | DEAP: Evolutionary Algorithms Made Easy[END_REF] is of interest for its improvement on the calculation of the crowding distance [START_REF] Fortin | Revisiting the NSGA-II crowding-distance computation[END_REF], and is therefore integrated in the framework.

Summary of the structure of the framework

How all of these components are then integrated in the optimization framework and communicate with each other can be described in Figure 1.

Each dash rectangle shows the limits of each component of the framework. The optimization process part contains all the evolutionary procedure of the GA. This optimization process sends an individual of the optimization to the evaluation function which use this individual to modify the DH network model (in our case modify the pipe diameters) and carry out a simulation. Simulation results are post-processed in the evaluation function to calculate the result of the objective function used in the optimization. The result of the objective function is then sent back to the optimization process and is affected to the individual. The evaluation function is a black box for the optimization process: it takes the individual as an argument and returns the value of the objective function.

Case study

In order to use the optimization framework on complex DH network sizing problems, we wanted to validate it on small-scale test cases first. This section describe the test case, which consist in a simple linear network and the formulation of the optimization problem.

DH network model description

The layout of the test case is displayed on Figure 2 the DH network consists of a straight horizontal two-tube network, with a set of N identical consumers equally distributed in the network. At the production level, we consider an ideal boiler along with a pump regulated to maintain a pressure head difference of 1 bar across each consumer of the network. The main model parameters are given in Table 1. Each consumer is modelled by its substation (SST). The substation model features an exchanger and a control valve regulated with a PID controller to provide the heat needed by the consumer. This toy problem is interesting to consider since, as we will show in section 4, it is possible to obtain an exact solution and to have a close understanding of the physics behind the results of the optimization. As the goal of this paper is to validate the optimization framework against an exact solution and to be able to keep an understanding of the physics, we opted for static simulation. In a preliminary step, we evaluated a naïve solution obtained by systematically considering a DN300 pipe for every sections of the toy problem. The main results of the simulation are provided in Table 2, and the solution is plotted in the next section's Figure 3. Those results show that those diameters are higher than needed to supply the heat into substation and are a consequent oversizing of the network.. 104 k€ For this test case, the optimization problem formulation is described in the following paragraphs.

Optimization formulation

The decision parameters are the pipe diameters of each section. Supply and return line of the same section of the network are supposed to be identical. Thus, the number of decision parameters is equal to the number of consumers, which is equal to 10.

The two objective functions for the optimization are the investment cost (CAPEX) needed for the pipes, and the operational pumping cost (OPEXpumping, as the part of OPEX related to pumping). The CAPEX is calculated by summing the cost for each pipe. The cost of each pipe € 𝑝𝑖𝑝𝑒 depends on its diameter and is calculated as € 𝑝𝑖𝑝𝑒 = 970 * 𝐷 + 30. According to DH Network operator, this linear relation is a fair approximation of piping prices.

As a widely accepted lifetime for the pipes is 20 years, we calculated the OPEX for 20 years of operation as well to make it comparable. According to French technical center CEREMA, DH networks operation during one year corresponds to 2500h at full power in average. Moreover, according to [START_REF] Li | District Heating Network Design and Configuration Optimization with Genetic Algorithm[END_REF], pumping cost can be approximated to 108€/MWh. Consequently, the OPEX on 20 years can be calculated with:

𝑂𝑃𝐸𝑋 20𝑦𝑟𝑠,𝑝𝑢𝑚𝑝𝑖𝑛𝑔 = ∑ 𝑂𝑃𝐸𝑋 𝑝𝑢𝑚𝑝𝑖𝑛𝑔 (1 + 0.05) 𝑦𝑟 20 𝑦𝑟=1
𝑤𝑖𝑡ℎ 𝑂𝑃𝐸𝑋 𝑝𝑢𝑚𝑝𝑖𝑛𝑔 = 𝑃 𝑝𝑢𝑚𝑝𝑖𝑛𝑔 * 2500 * 108 2 with 𝑃 𝑝𝑢𝑚𝑝𝑖𝑛𝑔 being the pumping power of the simulated network in MW. Please note that the choice for the values of pumping cost, actualization rate or hours of operation does not change which solution is optimal, but makes it comparable to the CAPEX as it reaches order of magnitude.

Constraints and parameters of the optimization

A single constraint is imposed on the DH network: at least 98% of the heat demand must be satisfied. If this constraint is not satisfied, the individual is eliminated from the optimization.

Parameters for NSGA2 algorithm are population size equal to 100 individual, 250 number of generation and a mutation rate of 0.1.

Decision space for the case study

The decision space consists of a selection of pre-insulated rigid steel pipes issued from a Logstor pipe catalogue. The pipes ranges from DN20 to DN350. Their rugosity height is 45µm; the insulation thermal conductivity is 0.0275 W/m.K. Other properties depends on the selected diameter and are detailed in Table 3.

Results

First, the process of building a reference Pareto Front used as an exact solution for this optimization problem is detailed. Then the results of the multi-objective optimization are presented and compared to the reference Pareto front.

Generation of a reference Pareto Front

We obtained the reference Pareto Front using a systematic but constrained search within the decision space. Every possible solution of the decision space respecting those constraints was calculated and only optimal solutions in the Pareto sense were kept. The constraint used to limit the computational burden was to only consider diameter sets which decrease with the distant from the supply unit. This simplification is sensible for this particular toy problem since the power to be transported by a section further from the boiler is lower than for the neighbour section closer to the boiler. The evaluation function part of the optimization framework was used to evaluate all of these sets. The reference Pareto front generated this way is plotted in red in Figure 3.

Results of the optimization

The output of the multi-objective optimization using the previously described framework and test case were obtained after a computation time of 16 hours on a regular laptop (Intel Core i7, 1.10GHz, with 8 cores allocated to the optimization). When analyzing the content of the solutions, all of them contain decreasing diameters while no constraints were imposed in the formulation of the optimization. This behavior confirm the relevance of the hypothesis used to generate the reference Pareto front. Diameters for the first pipe (close to the supply unit) range between DN200 and DN350 while diameters for the last pipe (further from the boiler and the pump) range between DN100 and DN200.

On the objective function values, the results of the optimization range from 193 k€ to 4 144 k€ for the pumping part of the OPEX calculated for 20 years of exploitation, and between 3 568 k€ to 6 424 k€ for the CAPEX.

All of the optimal solutions provided by the optimization are plotted in blue "x" on Figure 3 along with the reference Pareto front in red "+". A first qualitative comparison of the solutions provided in this figure is that the solutions seems to have a great diversity (the spread is similar when comparing obtained and reference Pareto set) and a very good convergence (the Pareto curves are close from one another).

The convergence of the optimization can be verified by the evolution of hypervolume indicator [START_REF] Zitzler | Performance assessment of multiobjective optimizers: An analysis and review[END_REF], which is available on Figure 4 and is compared on the same Figure to the hypervolume of the reference Pareto set. The stabilization of hypervolume at generation 110 implies that the optimization seems to reach a global optimum, as it does not discover optimum region of the objective space anymore. The global behavior of the search for optimal solution can be broken down into four phases, visible in the optimization thanks to Figure 4 and Figure 5. The first phase from the start to generation 70 is an approach of a local optimum. The second phase, until generation 120 is an increase of the diversity of the front with no major improvement in the convergence (there is still a small gap between the obtained and the reference Pareto front). Then until generation 180, the convergence is slightly improved and a few solutions are improved to the global optimum. No notable improvement are noticeable afterward.

The results of the optimization presented in this section, demonstrate a close correspondence to the reference results, in term of convergence to the same performance and in term of diversity of the solution proposed. Those results demonstrate the validity of the proposed method and framework to find accurate optimal solutions for the sizing of DH networks pipes on this test case and makes it possible to expand the use of the framework on further different problems.

Comparison to a local sizing method

In order to have a better understanding of the reference Pareto front, we propose to compare it to a solution elaborated with a local sizing method that is a solution widely employed in engineering offices. This solution can be computed by searching for the optimal diameter from the available diameters detailed in Table 3, for each section of the network. each section, the optimal diameter is the one that minimize the cost of the section. The global cost is determined by the aggregation the investment cost (CAPEX) and the pumping cost (OPEX) in a single value : we choose to calculate the Net Present Value (NPV) on 20 years. This economic indicator is widely used in the literature and is a key indicator for DH Network operators (Carpenè & Haeusler, 2019;[START_REF] Dufo-López | Multi-objective design of PV-wind-diesel-hydrogen-battery systems[END_REF][START_REF] Li | District Heating Network Design and Configuration Optimization with Genetic Algorithm[END_REF]. The NPV is computed as follows:

𝑁𝑃𝑉 = 𝐶𝐴𝑃𝐸𝑋 + ∑ 𝑂𝑃𝐸𝑋 (1 + 0.05) 𝑦𝑟 20 𝑦𝑟=1 4
The OPEX and the CAPEX are computed with the same procedure as the one detailed in section 3 for the multiobjective optimization..The pumping power of the OPEX is evaluated by calculating the regular pressure loss across the section. This term is computed with the same correlation as in DistrictLab-H, ie. proposed by [START_REF] Swamee | Explicit Equations for Pipe-Flow Problems[END_REF]. In this approach, the problem is simplified to a set of mono-objective, mono-variable minimization problems.

The optimal solution of the local sizing method is compared to the closest solution from the reference Pareto set. This solution is selected by performing the same aggregation of the OPEX and CAPEX in order to pick the optimal solution of the same mono-objective problem. The comparison of the diameters provided by both of those methods is available in Table 4. 125 Proposed diameters are the same for the section 1 to 6, and after those sections, the local sizing method result in diameters smaller than the one present in the reference Pareto set. This suggest that in the optimization framework, the OPEX has more influence on the diameters selected than in the local sizing method, and thus that the pumping power is under estimated in the local sizing method.

DistrictLab-H simulation's results for regular pressure loss in the pipes are in accordance with the values proposed by the local sizing method, but other pressure losses are simulated by DistrictLab-H in other components such as substations and control valve. These components have an effect on the system global operating cost that is not taken into account when sizing each section independently. Indeed, system level effects arise due to dependency of the pumping power to head losses with respect to the global mass flow rate, whereas local sizing only accounts for the mass flow rate effectively circulating in the pipes of interest. Thus we think that the hypothesis of independence of each section of the DH network does not enable an optimal sizing of the pipes.

Even if further research is needed to understand fully the difference between the solutions given in the local sizing method and the optimization framework, those differences suggest that a DH network is a complex system that can benefit from a sizing method as the one proposed by our framework that treats the network as a whole and which accounts for inner interactions.

Conclusion

This paper presented a sizing framework for distribution pipes in DH networks. The framework is based on multi-objective optimization and its components were presented in section 2. A toy problem was used to test and validate the framework on a small scale problem. The results section of this paper suggests that the optimization framework is capable of closely approximating enough the Pareto front of an optimization for our application. It can thus be used for decision support regarding optimal sizing for a DH network distribution, at least for simple network layout. First, the validation against a reference Pareto front enables us to ensure the quality of provided results on a simple test case and in the capacity of the framework to tackle this type of problems. Secondly, the comparison to a local sizing method show an interest to the method of sizing integrated in our optimization framework.

The choices made in the development of this framework lean towards the resolution of more complex problems: the introduction shows us that problems with more than a hundred parameters were solved with similar methodological choices [START_REF] Hohmann | Multi-objective optimization of the design and operation of an energy hub for the EMPA campus[END_REF][START_REF] Li | District Heating Network Design and Configuration Optimization with Genetic Algorithm[END_REF][START_REF] Vesterlund | Design Optimization of a District Heating Network Expansion, a Case Study for the Town of Kiruna[END_REF]. With the small-scale validation, we are confident to tackle such problems, and to use the DistrictLab-H simulation platform for dynamic simulation in these optimization problems while keeping a reasonable calculation cost. Dynamic simulation will make it possible to optimize the sizing of the network while taking into account the possibilities of storage, to include the impact of transient states on the operation cost of the network.

Main limitations of this optimization framework lie in two points. First, the complexity of the evaluation function implies that the real Pareto front is often too costly to find, and that the Pareto front obtained as a results of the optimization is an approximated Pareto Front. Therefore, solutions of the obtained Pareto Front may be used as a good starting point for further studies before the implementation of those solutions on site. The second limitation lies into the assessment of the quality of the results of multi-objective optimization. This issue is tied with the calculation cost of such an optimization: the choice of the number of evaluation has to be made carefully to ensure the convergence while limiting the number of evaluation. Further work will also address the validation of the framework on other toy problems, such as a tree-shaped DH network and a loop shaped DH network before using the framework on a city-scale complex problem.

Figure 1 :

 1 Figure 1 : Structure and components of the optimization framework

Figure 2 :

 2 Figure 2 : Layout of test case 1

Figure 3 :

 3 Figure 3 : Pareto fronts for the linear test case : comparison between the results of the optimization framework and the reference In order to have a comprehensive metric to compare both of the Pareto fronts, we chose to calculate an indicator representing the average distance in euros between each point of the obtained Pareto Front and its nearest neighbour in the reference Pareto front. With 𝐴 = {𝑎 1 , 𝑎 2 , … , 𝑎 |𝐴| } being the obtained Pareto front and its objective function value, and 𝑍 = {𝑧 1 , 𝑧 2 , … , 𝑧 |𝑍| } being the reference Pareto front and its objective function value, the indicator is computed like the generational distance (Van Veldhuizen, 1999):

Figure 4 :

 4 Figure 4 : Evolution of the hypervolume of the obtained Pareto front along the optimization

Figure 5 :

 5 Figure 5 : Evolution of generational distance between reference and obtained Pareto front along generations, logarithmic scale

Table 1 :

 1 model parameters for test case with 10 consumers

	Sizing total heat demand	50 MW
	Total length of DH network	10km supply line + 10 km return line
	Boundary condition	P=2 bars on return line at boiler
	Temperature setpoint of the boiler	135°C

Table 2 :

 2 main simulation results for the "naïve" solution with DN300 pipes everywhere

	Total heat demand in substations	50 MW
	Total heat injected by heat generator	50,8 MW
	Total heat loss	700 kW
	Total pumping power (𝑃 𝑝𝑢𝑚𝑝𝑖𝑛𝑔 )	69 kW
	CAPEX	6 420 k€
	OPEX20yrs,pumping

Table 3 :

 3 Decision space : possible pipe specifications

	DN	20	25	32	40	50	65	80
	Insulation thickness (mm) 31	28	34	31	32	32	36
	DN	100	125	150	200	250	300	350
	Insulation thickness (mm) 42	42	41	48	63	63	72

Table 4 :

 4 comparison of diameters for the optimal solution given by the local sizing method and the optimization framework

	Section	1	2	3	4	5	6	7	8	9	10
	Local sizing	250	250	250	250	200	200	150	150	125	80
	Ref Pareto set 250	250	250	250	200	200	200	200	150

Corresponding author. Email address : roland.baviere@cea.fr


Jean-François Lemaître
email: jean-francois.lemaitre@univ-lyon1.fr

Louise Lemaître

Frédéric Cheynel

G Douhard

François Bourgoin

Hubert Débias

Emmanuelle Ferté

Sylvia Gilot-Fromont

Maryline Pardonnet

Benjamin Pellerin

Rey

Louise Cheynel

| Frédéric Douhard

| Gilles Bourgoin

François Débias

| Hubert Ferté

Emmanuelle Gilot-Fromont

Sylvia Pardonnet

| Maryline Pellerin

Benjamin Rey

| Cécile Vanpé

| A J Mark Hewison

Jean-Michel Gaillard

The influence of early-life allocation to antlers on male performance during adulthood: Evidence from contrasted populations of a large herbivore

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.

2. According to the principle of allocation, high energy allocation to growth or sexual traits in males should lead to reduced energy allocation to the maintenance of cellular and physiological functions, potentially causing them to age faster, with impaired survival.

3. We evaluated the short-term and delayed consequences of energy allocation to antlers early in life in two contrasted populations of roe deer, Capreolus capreolus.

Although most males mate successfully for the first time in their fourth year, antlers are grown annually from the first year of life onwards. We tested the prediction that a high level of allocation to antler growth during the first two years of life should lead to lower body mass, antler size and survival during the early and late prime stages, as well as to reduced longevity overall.

4. Growing and carrying long antlers during the first years of life was not associated with any detectable cost in the late prime stage. The positive association between antler growth in early life and adult body mass instead supports that fawn antler acts as an honest signal of phenotypic quality in roe deer. 5. For a given body mass, yearling males growing longer antlers displayed impaired performance during their late prime. We also found a trend for a short-term survival cost of allocation to relative antler length during the second year of life.

| INTRODUCTION

The evolution of complex secondary sexual traits in males is generally explained by the positive relationship between the intensity of their expression and breeding success [START_REF] Andersson | Sexual selection[END_REF], as is the case for plumage coloration in many birds [START_REF] Brommer | The colour of fitness: Plumage coloration and lifetime reproductive success in the tawny owl[END_REF]. All models of sexual selection that have been developed thus far [START_REF] Kuijper | A guide to sexual selection theory[END_REF] share the common and mandatory assumption that growing and maintaining secondary sexual traits is costly [START_REF] Kotiaho | Costs of sexual traits: A mismatch between theoretical considerations and empirical evidence[END_REF], although these costs may vary in relation to individual condition [START_REF] Rowe | The lek paradox and the capture of genetic variance by condition dependent traits[END_REF]. However, measuring the fitness costs of bearing extravagant sexual traits empirically is challenging and most studies of life-history (e.g. higher mortality due to higher predation risk) or physiological (e.g. lower performance of the immune system) costs to date have reported contrasting patterns [START_REF] Kotiaho | Costs of sexual traits: A mismatch between theoretical considerations and empirical evidence[END_REF]. A meta-analysis across birds, spiders, insects and fish revealed that the expression of secondary sexual traits is strongly condition-dependent as the expression of a given sexual trait (known to predict mating success) is generally positively associated with longevity [START_REF] Jennions | Sexually selected traits and adult survival: A meta-analysis[END_REF]. However, while males in good condition during the early stages of life are able to develop the most conspicuous sexual traits, they may suffer from a higher degree of somatic deterioration compared to males that are in low condition during early life, which can lead to an accelerated rate of ageing [START_REF] Adler | Phenotypes optimized for early-life reproduction exhibit faster somatic deterioration with age, revealing a latent cost of high condition[END_REF].

In line with these results, it has recently been argued that the costly nature of secondary sexual traits should be rooted in the theories currently proposed to explain the evolution of ageing [START_REF] Hooper | Sexual competition and the evolution of condition-dependent ageing[END_REF][START_REF] Tidière | Does sexual selection shape sex differences in longevity and senescence patterns across vertebrates? A review and new insights from captive ruminants[END_REF]. O ur under s t anding of the ageing process t akes it s root in the principle of allocation [START_REF] Cody | A general theory of clutch size[END_REF], a cornerstone of life-history theory [START_REF] Stearns | The evolution of life histories[END_REF]. This principle states that individuals need to share a finite pool of resources among competing functions such as growth, reproduction and survival. When it concerns a long-lasting or delayed cost of a given allocation, this principle is closely related to the disposable soma theory of ageing [START_REF] Kirkwood | The evolution of aging and longevity[END_REF]. This theory postulates that individuals that allocate a substantial amount of resources to reproduction will have less resources to devote to the set of molecular mechanisms (e.g. DNA repair mechanisms) involved in somatic maintenance [START_REF] Kirkwood | Evolution of senescence: Late survival sacrificed for reproduction[END_REF]. Therefore, individuals that allocate much to reproduction early in life will accumulate a variety of damages that will ultimately cause a shortened life span or an increased rate of reproductive senescence (Lemaître & Gaillard, 2017). The central role played by early-vs. late-life trade-offs in the evolution of ageing constitutes the core of the life-history theory of ageing [START_REF] Kirkwood | The Disposable soma theory: Origins and evolution[END_REF] and has been supported by many empirical studies [START_REF] Hammers | The impact of reproductive investment and early-life environmental conditions on senescence: Support for the disposable soma hypothesis[END_REF][START_REF] Jankowiak | Age-specific reproduction and disposable soma in an urban population of Common Blackbirds Turdus merula[END_REF]. However, most of these studies have focused on females [START_REF] Lemaître | Early-late life trade-offs and the evolution of ageing in the wild[END_REF] and there is no consensus yet on whether early-vs. late-life trade-offs are prevalent in males.

During the reproductive season, males devote a substantial part of their resources to sexual competition [START_REF] Andersson | Sexual selection[END_REF]. While sexual competition can take different forms, male breeding success is determined by allocation to secondary sexual traits in most species [START_REF] Bonduriansky | Sexual selection and allometry: A critical reappraisal of the evidence and ideas[END_REF][START_REF] Emlen | The evolution of animal weapons[END_REF].

Several studies have revealed that the allocation to secondary sexual traits can be associated with diverse physiological costs [START_REF] Garratt | Oxidative stress and conditiondependent sexual signals: More than just seeing red[END_REF]. In the Australian painted dragon (Ctenophorus pictus), males that maintained a high level of head coloration suffered from a rapid rate of telomere loss [START_REF] Giraudeau | Ageing and the cost of maintaining coloration in the Australian painted dragon[END_REF]. This finding suggests that males with particularly conspicuous traits might have impaired fitness later in life. While current theories predict that males should suffer from long-term or delayed costs of growing costly sexual traits, such allocation can levy immediate fitness costs. Thus, successful reproduction often has shortterm reproduction and survival costs in male birds [START_REF] Bleu | Reproductive costs in terrestrial vertebrates: Insights from bird studies[END_REF], which might potentially be driven by strong allocation to sexual traits.

Overall, measuring costs of allocation to secondary sexual traits during early life is important because the expression of secondary sexual traits varies over the lifetime (Lemaître & Gaillard, 2017), and also because allocation to secondary sexual traits early in life increases the absolute level of allocation to the entire growth process, which potentially leads to detrimental fitness consequences later in life [START_REF] Lee | Experimental demonstration of the growth rate -lifespan trade-off[END_REF][START_REF] Metcalfe | Growth versus lifespan: Perspectives from evolutionary ecology[END_REF]. This avenue of research should boost our understanding of the relationships linking individual development with population dynamics. So far, the few studies that assessed the long-term fitness consequences of early allocation to secondary sexual traits in mammalian males have mostly been performed on bovids (see Table 1 for a review), and showed that the intensity of early horn growth is rarely associated with a delayed survival cost (Table 1). While the horns of bovids grow continuously over most of their lifetime, the amount of energy allocated to the annual horn growth is quite limited. On the contrary, male cervids grow antlers that are cast 6. Yearling males that grow long antlers relative to their mass might display a fast lifehistory tactic. We argue that differential allocation to secondary sexual traits generates a diversity of individual trajectories that should impact population dynamics.

KEYWORDS

ageing, Capreolus capreolus, life history, roe deer, secondary sexual traits, senescence and fully regrown every year [START_REF] Chapman | Antler bones of contention[END_REF]. This repeatedly high annual allocation of resources to antler growth is thus expected to be particularly costly. So far, the only study that quantified delayed costs of growing antlers reported that male red deer (Cervus elaphus) that allocated substantially to antler growth between 4 and 9 years of age did not suffer any detrimental consequences on harem size or antler length later on [START_REF] Lemaître | Early life expenditure in sexual competition is associated with increased reproductive senescence in male red deer[END_REF]. However, young males grow antlers from their first year of life, well before they can successfully mate and at a time when they must allocate intensively to body growth [START_REF] Douhard | The cost of growing large: Sex-specific costs of post-weaning growth on body mass senescence in a wild mammal[END_REF]. Furthermore, delayed fitness costs of early allocation to antler growth might impact survival at the end of the prime-aged adulthood (i.e. when the onset of actuarial senescence occurs; [START_REF] Gaillard | Senescence in mammalian life-history traits[END_REF] rather than reproductive success (Lemaître & Gaillard, 2017). Surprisingly, the long-term fitness consequences of early development of secondary sexual traits have not yet been investigated. Moreover, while there is increasing evidence that early-life environmental conditions mediate trade-offs between early-and late-life performance [START_REF] Balbontín | Environmental conditions during early-life accelerate the rate of senescence in a short-lived passerine bird[END_REF], whether the interplay between environmental conditions and the expression of secondary sexual traits early in life shapes the magnitude of delayed fitness costs remains unknown.

We propose to fill this knowledge gap by investigating the consequences of allocation to antlers in early life (defined here as the first two years of life) on survival, body mass and antler size during adulthood in two populations of European roe deer (Capreolus capreolus) living in contrasted habitats. In this slightly polygynous mammal, males are territorial from early spring until mid-July-late August and generally begin mating during their fourth year, although some males can sire offspring in their third year of life [START_REF] Vanpé | Age-specific variation in male breeding success of a territorial ungulate species, the European roe deer[END_REF].

Male roe deer systematically grow antlers each year from their first year of life onwards. Growing large antlers in early life well before mating opportunities occur should thus come at a cost in terms of survival, body mass or antler growth later in life. Finally, as resource availability mediates the trade-off between growth and somatic maintenance, such costs should be more pronounced in populations experiencing poor environmental conditions.

| MATERIALS AND METHODS

| Study population

We studied the populations of roe deer at Trois-Fontaines (TF) and summers. This site has rich soils and provides high-quality habitat for roe deer (Pettorelli et al., 2006). In contrast, in Chizé (2,614 ha), located in western France (46°50′N, 0°25′W), the climate is temperate oceanic with Mediterranean influences. This forest has low productivity due to poor-quality soils and frequent summer droughts (Pettorelli et al., 2006) and thereby provides a quite poor habitat for roe deer in most years.

Roe deer have been intensively monitored using a long-term capture-mark-recapture programme since 1975 and 1977 for Trois-Fontaines and Chizé, respectively. Every year and for each site, 10-12 days of capture is organized between December and March (Gaillard et al., 1993a). Once an individual is captured, its sex and body mass (to the nearest 50 g) are recorded and a basic clinical examination is performed.

| Allocation to antler growth during the first two years of life

Antlers grow each year from the end of November until early March (i.e. before the velvet is shed). Roe deer are in hard antler until the end of October/early November when they cast their antlers [START_REF] Sempéré | Relationship between antler development and plasma androgen concentrations in adult roe deer[END_REF]. During their first year of life, males develop a small "button" (i.e. a very short antler with a round top instead of a spike) on top of their pedicle before growing their first full antlers at 1 year of age (Figure 1). Then, males generally develop a four-point (in yearlings) or five-/six-point (in adults) head. Antler length was measured to the nearest 0.5 cm along the external side of the main beam, from the base of the antler to the top of the main beam [START_REF] Vanpé | Antler size provides an honest signal of male phenotypic quality in roe deer[END_REF]. Antler length was used as a proxy for allocation to antler growth because there is a strong covariation among the different antler traits in deer [START_REF] Lemaître | Early life expenditure in sexual competition is associated with increased reproductive senescence in male red deer[END_REF]. When an individual was caught more than once in a given year, only the last measure was considered. When both antlers were intact, the average length was computed and retained for the analyses. When one antler was broken, the value of the remaining intact antler was retained. When both antlers were broken, the individual was removed from the dataset. Although a few captures occurred in October-November in the very first years of monitoring, we restricted our analyses to data collected from 1 December to avoid possible overlap between antlers growing during the current year and uncast antlers from the previous year. All males included in our analyses were of known age because they were caught within their first year of life, when age can be assigned without error from tooth eruption and wear patterns [START_REF] Hewison | Tests of estimation of age from tooth wear on roe deer of known age: Variation within and among populations[END_REF].

In roe deer, the antler growth cycle is governed by complex interactions between the photoperiod and the endogenous rhythm of the hypothalamic-pituitary-gonadal axis [START_REF] Sempéré | Temporal and seasonal relationships between LH, testosterone and antlers in fawn and adult male roe deer[END_REF].

Therefore, for a given cohort, the date of capture provides a better predictor of antler size during antler growth than birth date [START_REF] Vanpé | Antler size provides an honest signal of male phenotypic quality in roe deer[END_REF]. We analysed the relationship between antler length and the Julian date of capture (with 1 December as day 0). As antler growth is a nonlinear process [START_REF] Goss | Deer antlers: Regeneration, function and evolution[END_REF] then increasing further, but with a different slope). To estimate these threshold values, the deviance profile of the models including the slopes was used and the date of capture providing the lowest deviance was selected as the threshold value [START_REF] Ulm | On the estimation of the threshold values[END_REF]. Model selection was based on the Akaike Information Criterion (AIC). We calculated AIC weights (AICw) to assess the relative likelihood that a given model was the best among all the fitted models [START_REF] Burnham | Model Selection and multimodel inference: A practical information-theoretic approach[END_REF]. We selected the model with the lowest AIC except when the difference in AIC (ΔAIC) between two competing models was <2, in which case we retained the simplest model. Once the model that best captured the shape of the relationship was identified, we looked for potential differences between populations by including "population" as a fixed factor, either as an additive effect or in interaction with the Julian date.

As previously observed [START_REF] Vanpé | Antler size provides an honest signal of male phenotypic quality in roe deer[END_REF], the linear model best represented variation in antler length (log-transformed) with date of capture in fawns (Table S1; Figure 2a; slope of 0.022 ± 0.002, N = 224). Fawn antler length did not differ between populations (∆AIC of + 1.93 for the model including "population"), and the daily growth rate did not differ between populations (∆AIC of +3.57 for the model including the "population × date of capture" interaction term). We computed a standardized fawn antler length from the selected model by adjusting antler length to the median date of capture (i.e. 8 February) and pooling the data for both populations.

The best model accounting for variation in yearling antler length in relation to capture date was the threshold model with a single slope (Table S2). Antler length increased linearly with date of capture until day 59 (i.e. 28 January) (Figure 2b; slope of 0.024 ± 0.002, N = 257) and remained constant thereafter. Antler length did not differ between populations (∆AIC of +1.81 for the model including "population"), but the threshold date did (∆AIC of -5.17 for the model including the "population x date of capture" interaction term). In both populations, antler length increased linearly with the date of capture (0.032 ± 0.004, N = 100, and 0.021 ± 0.002, N = 157 at Trois-Fontaines and Chizé, respectively), but the threshold date was earlier at Chizé (day 46 corresponding to 16 January; term effect of allocation to fawn relative antler length (see statistical analysis), we thus standardized body mass by the Julian date (8 February), as done previously [START_REF] Vanpé | Antler size provides an honest signal of male phenotypic quality in roe deer[END_REF]. We followed the same procedure for yearling antler length. As expected, there was no detectable variation in body mass in relation to date at this age (Figure S2; slope of 0.0003 ± 0.0005, N = 257), and the yearling body mass was thus not corrected for date of capture.

| Covariation between fawn and yearling antler length

We used two metrics to describe the level of allocation to antler growth during early life: fawn and yearling absolute standardized antler length. As we expected heavier males to allocate more resources to antler growth in absolute terms than lighter males [START_REF] Vanpé | Antler size provides an honest signal of male phenotypic quality in roe deer[END_REF], we measured the effect of relative allocation to antlers during the first two years of life by adding body mass (i.e. fawn body mass standardized for the date of capture or yearling non-standardized body mass) as a covariate in the corresponding set of models (see below).

Yearling antler size could potentially be influenced by the allocation to antlers the year before, and we thus first tested whether fawn and yearling antler length could be considered as independent.

As the physiological pathway governing antler growth is the same across all ages [START_REF] Sempéré | Temporal and seasonal relationships between LH, testosterone and antlers in fawn and adult male roe deer[END_REF], we expected a positive relationship to occur between absolute antler length during the first and second years of life, and we thus tested for a positive correlation between allocation to fawn and yearling antler growth using a linear regression. We replicated the analyses on relative (to body mass) measures of antler length. To do this, we did not include fawn and yearling body mass (log-transformed) to avoid multicollinearity issues. We thus first computed fawn relative antler length as the residuals from the linear regression between absolute standardized antler length and fawn standardized body mass (slope of 1.98 ± 0.37, N = 224; Figure 3a). We followed the same procedure to compute yearling relative antler length (slope of 0.68 ± 0.17, N = 257; Figure 3b). We then computed yearling relative antler length as the residuals from this linear regression between yearling relative antler length and body mass. Both fawn and yearling antler length were measured on 93 males (76 in Chizé and 17 in Trois-Fontaines).

Measures of fawn and yearling antler length were independent (absolute metrics: slope of 0.043 ± 0.033; Figure 3c, Table S3; relative metrics: slope of 0.001 ± 0.04; Table S4).

| Early and late prime-age adulthood traits

To assess the consequences of male allocation to antler growth during early life on performance during early and late prime-age adulthood, we used three survival metrics: the probability of reaching 3 years of age (survival metric for short-term survival), the probability of reaching 6 years of age (survival metric for long-term survival) and longevity. We also used two phenotypic traits that are positively associated with fitness in roe deer: adult body mass [START_REF] Gaillard | Body mass and individual fitness in female ungulates: Bigger is not always better[END_REF] and adult antler length [START_REF] Vanpé | Assessing the intensity of sexual selection on male body mass and antler length in roe deer Capreolus capreolus: Is bigger better in a weakly dimorphic species?[END_REF]. For a given male, adult body mass and antler length were measured as the median values of measures recorded between 4 and 6 years of age (body mass) and 4 and 8 years of age (antler length), as these traits subsequently show senescence [START_REF] Douhard | The cost of growing large: Sex-specific costs of post-weaning growth on body mass senescence in a wild mammal[END_REF][START_REF] Vanpé | Antler size provides an honest signal of male phenotypic quality in roe deer[END_REF]. Longevity (from 1.5 to 14 years) was measured as the exact number of years an individual lived, including only those individuals that were intensively monitored by observations and that died from natural causes (i.e. individuals that were removed from the population for translocation or that died from human-related incidents such as car collisions and capture-related mortality were right-censored).

| Analysis of the consequences of strong allocation to antlers during early life

We investigated potential trade-offs between allocation to antlers during early life and both condition and survival during adulthood and possible between-population differences in these trade-offs.

For each combination of late-prime trait (dependent variable) and absolute antler length in early life, a set of competing models was fitted (e.g. Table S5a). When testing the possible long-term deleterious consequences of allocation to relative antler length early in life, body mass (i.e. fawn standardized body mass or yearling nonstandardized body mass) was added as a covariate in all models of the set. All models included the cohort (year of birth) fitted as a random effect (32 cohorts between 1975 and 2006). When analysing adult antler length, we also included adult body mass as a covariate (e.g. Table S5b) to measure the consequences of early allocation to antler growth on both absolute and relative adult antler length. For analysis of longevity, adult body mass and adult antler length, we fitted linear mixed-effects models ("lmer" function in the r package lme4, [START_REF] Bates | Fitting linear mixed-effects models using lme4[END_REF]. For survival beyond 3 and 6 years of age, we used generalized linear mixed-effects models with survival entered as a binomial variable. We tested the short-term survival costs (i.e. survival beyond 3 years of age) only for yearling antler length as our sample size was too limited to investigate short-term survival costs in fawns.

All traits except survival were log-transformed. For each set of analyses, we used AIC for model selection (see above for a full description of the procedure). All analyses were performed using r version 3.4.0 [START_REF] Core | R: A language and environment for statistical computing[END_REF], and all results are reported as M ± SE.

In our dataset, two individuals had particularly short antlers for their age (one fawn with antlers measuring 1 mm and one yearling with antlers measuring 10 mm). All analyses were performed with and without these males, and in most cases, results were qualitatively similar (Table S5 vs. Table S6 for fawn antler length;

Table S7 vs. Table S8 for yearling antler length). As we did not have any biological reason to exclude these males, we decided to present the results obtained with the complete dataset. We explicitly specified in the text when including these outliers mattered. Population-specific sample sizes used for each analysis are provided in Table 2.

| RESULTS

| Late prime-age adulthood consequences of allocation to fawn antler growth

Absolute and relative antler length were not associated with longevity or survival beyond 6 years of age (longevity: slope of 0.04 ± 0.05 and slope of 0.05 ± 0.06; survival: slope of 0.07 ± 0.20 and slope of -0.71 ± 0.47 for absolute and relative antler length, respectively).

In all cases, the constant model was selected (Table S5a andS5b). S5a, S5b), the selected model only included a positive effect of body mass (slope of 0.39 ± 0.06; Figure S3). The selected model for adult antler length only included a positive effect of adult body mass (slope of 0.61 ± 0.22, N = 111), as expected for this allometric relationship (Figure S4).

| Early and late prime-age adulthood consequences of yearling allocation to antler growth

No relationship occurred between allocation to absolute antler length and the probability of reaching 3 years of age (slope of -0.56 ± 0.48; Table S7a). On the other hand, yearling tended to suffer from carrying long antlers for their body mass (slope of -0.90 ± 0.55), even if the most parsimonious model included body mass only (Table S7b).

Heavy individuals survived better (slope of 2.63 ± 1.29). Increasing by 1 kg the average body mass was associated with an increased survival of 2.5%.

Allocation to absolute antler length was not associated with any longevity or long-term survival cost in either roe deer population (slopes of -0.13 ± 0.11 and of -0.29 ± 0.35 for longevity and survival beyond six years, respectively; see Table S7a andS7b).

Analyses of relative antler length suggested that antler length for a given mass was negatively associated with longevity (Table 3, Figure 5). Increasing the average antler size by 5 cm at a given mass led to a reduction in longevity by 4.84%. However, most models performed equally well (Table S7b). For survival beyond six years, S7b). At Chizé, a strong relative allocation to antler growth led to a decrease in the probability of surviving beyond six years of age (Table 4, Figure 6), while relationship occurred at Trois-Fontaines (Table 4, Figure 5c,d). When the individual with 10-mm antlers was removed from the analysis, a negative effect of relative antler length was detected across both populations (Table S8b). Increasing the average relative antler size by 1 cm led survival to be reduced by 2.4%. Both absolute and relative antler length were independent of adult body mass (slope of 0.03 ± 0.02 and slope of -0.02 ± 0.01 for absolute and relative antler length, respectively; Table S7a andS7b). The selected model for absolute antler length only included population. Males were heavier at Trois-Fontaines than at Chizé, independently of antler growth. The selected model for relative antler revealed that body mass was positively associated with adult body mass (slope of 0.62 ± 0.05, N = 126). Finally, for adult antler length, the selected model contained both standardized antler length and body mass, revealing that relative antler length is positively associated with long antlers during adulthood (slope of 0.18 ± 0.08; Table S7a andS7b).

| DISCUSSION

Contrary to our predictions, we did not detect any evidence of delayed costs of allocation to antler growth during the first year of life in two roe deer populations living in highly different ecological contexts. Fawns with the longest (absolute or relative) antlers also had the highest body mass as adults, suggesting that the length of the first antlers reliably indicates individual performance during adulthood. This interpretation is supported by our finding that fawn antler length (but also yearling) was also positively associated with body mass, which indicates that antler length is a condition-dependent trait [START_REF] Vanpé | Antler size provides an honest signal of male phenotypic quality in roe deer[END_REF]. In yearling, the growth of long absolute antlers has no influence on future male performance. However, for Although the exact physiological mechanisms regulating the full antler growth cycle are still debated [START_REF] Bartos | Endocrine relationships between rank-related behaviour and antler growth in deer[END_REF][START_REF] Price | Exploring the mechanisms regulating regeneration of deer antlers[END_REF], experimental manipulation of antler size has revealed that the annual growth and casting of antlers is governed by a complex interaction between photoperiod and several hormones, especially insulin-like growth factor 1 (IGF-1) and testosterone [START_REF] Price | Exploring the mechanisms regulating regeneration of deer antlers[END_REF]. While the relative contribution of these hormones to antler size is difficult to decipher [START_REF] Ditchkoff | Concentrations of insulin-like growth factor-I in adult male whitetailed deer (Odocoileus virginianus): Associations with serum testosterone, morphometrics and age during and after the breeding season[END_REF], positive correlations between levels of testosterone and IGF-1 and antler length have been repeatedly reported [START_REF] Bartos | Endocrine relationships between rank-related behaviour and antler growth in deer[END_REF], notably in roe deer [START_REF] Schams | Antler growth in roe deer: Effects of casting and correlations with growth hormone and insulin like growth factor I[END_REF][START_REF] Sempéré | Relationship between antler development and plasma androgen concentrations in adult roe deer[END_REF]). An increasing number of studies have highlighted the potential deleterious effects of these hormones [START_REF] Foo | The effects of sex hormones on immune function: A meta-analysis[END_REF], which are assumed to play a pivotal role in governing life-history trade-offs [START_REF] Dantzer | Mediation of vertebrate life histories via insulin-like growth factor-1[END_REF], notably between allocation to sexual competition during early life and ageing [START_REF] Brooks | Life history evolution, reproduction, and the origins of sex-dependent aging and longevity[END_REF]. As yearlings with the longest relative antler length also carry long antlers as adults, repeated exposure to high hormonal levels might lead to reduced survival. The relationship between yearling and adult antler size suggests that yearlings with disproportionally large antlers will also have high competitive abilities over their entire life. Roe deer males are strongly territorial from early spring (March-April) to late August-early September [START_REF] Andersen | The European roe deer: The biology of success[END_REF], so that the costs of territory defence are likely higher than the costs associated with antler growth. Our results suggest that these costs might be higher when environmental conditions are poor, al- reproduction is positively associated with longevity in mammalian females (e.g. [START_REF] Descamps | Best squirrels trade a long life for an early reproduction[END_REF], and our results suggest that this might also be the case in males.

Yearlings that also grow relatively long antlers might suffer from more harassment by territorial males that are already established than those with smaller weapons [START_REF] Wahlström | The significance of male-male aggression for yearling dispersal in roe deer (Capreolus capreolus)[END_REF], as suggested by the detected tendency for short-term survival costs.

A few recent studies have revealed that high levels of sexual competition can have long-term consequences in terms of age-specific decline in individual performance [START_REF] Beirne | Sex differences in senescence: The role of intra-sexual competition in early adulthood[END_REF][START_REF] Lemaître | Early life expenditure in sexual competition is associated with increased reproductive senescence in male red deer[END_REF]. More generally, one promising avenue of research linking functional traits and population performance would be to investigate the dynamics of antler length very early in life and its impact population growth. Our study highlights that growing large secondary sexual traits early in life might be a reliable indicator of individual performance, suggesting that variation in the expression of sexual traits within populations might directly influence population demography. Integral projection models (IPMs) provide a powerful tool to investigate this question [START_REF] Coulson | Integral projections models, their construction and use in posing hypotheses in ecology[END_REF], but so far, IPMs have mostly focused on body mass or size as the focal trait [START_REF] Vindenes | Individual heterogeneity in life histories and eco-evolutionary dynamics[END_REF] and most applications have been restricted to females. The recent development of two-sex models [START_REF] Schindler | Sex-specific demography and generalization of the Trivers-Willard theory[END_REF] now offers a way to quantify how variation in the expression of sexual traits such as antler length influences the demography of a population, provided the required data (i.e. sex-specific relationships between demographic parameters [recruitment and survival] and antler length, father-offspring relationships and age-specific changes of antler length, pedigrees) will be available. As paternity analyses are becoming very common in wild vertebrates, we hope that our study will stimulate such studies that might reveal important insights on the role of sexual selection on population demography.

Finally, one striking paradox from our study is that the honest signal of individual performance conveyed by antler length is apparent in fawns when males do not compete for females. A possible explanation for the growth of long antlers in the pre-reproductive stage is the hormonal regulation of antler development that might be linked to physiological pathways involved in health and somatic maintenance. These constraints may have prevented natural selection from counterselecting the growth of substantial antlers in cervids well before they serve their purpose in male-male sexual contests [START_REF] Sempéré | Temporal and seasonal relationships between LH, testosterone and antlers in fawn and adult male roe deer[END_REF]. Future work should test whether long antlers in the first years of life reflect the quality of the ejaculate and thus, possibly, male breeding success [START_REF] Malo | Antlers honestly advertise sperm production and quality[END_REF]. Whether this signal is already present in fawn antlers is unknown (but note that young roe deer with long antlers also have large testes; [START_REF] Wahlström | The significance of male-male aggression for yearling dispersal in roe deer (Capreolus capreolus)[END_REF]. This is the missing piece of the jigsaw for our understanding of the relationships between age-specific functional traits, health and demographic parameters.

  Chizé (CH) living in enclosed forests. Trois-Fontaines (1,360 ha), located in north-eastern France (48°43′N, 4°55′E), has a continental climate characterized by moderately severe winters and warm rainy

FIGURE 1 A

 1 FIGURE 1 A fawn with visible buttons (photograph reproduced with permission from F. Débias)

Figure 2d )

 2d Figure 2d) than at Trois-Fontaines (day 72 corresponding to 11February; Figure2c), which matches the longer day length at Chizé than at Trois-Fontaines in January-February. We then computed a standardized yearling antler length by adjusting antler length to the median date of capture (i.e. 8 February) using the model selected for each population.

FIGURE 2

 2 FIGURE 2 Relationship between antler length (log-transformed) during the first two years of life and Julian date of capture. Fawn antler length (at 8 months of age) for males at Chizé and Trois-Fontaines (a); yearling antler length (at 20 months of age) at Chizé and Trois-Fontaines (b); yearling antler length (at 20 months of age) at Trois-Fontaines (c); and yearling antler length (at 20 months of age) at Chizé (d)

FIGURE 3

 3 FIGURE 3 Relationship between fawn standardized antler length (log-transformed, measured at 8 months of age) and fawn standardized body mass (log-transformed) (a). Relationship between yearling standardized antler length (log-transformed, measured at 20 months of age) and yearling body mass (log-transformed, measured at 20 months of age) (b). Relationship between yearling standardized antler length (log-transformed and measured at 20 months of age) and fawn standardized antler (log-transformed, measured at 8 months of age)

  However, when analysing the effect of absolute antler length on longevity, the selected model contained population as a fixed factor. For a median value of antler length, males lived about one year longer at Chizé than at Trois-Fontaines (Trois-Fontaines: 4.49 years 95% CI [3.78-5.33], N = 43; Chizé: 5.39 years 95% CI [4.80-6.05], N = 94, estimates averaging all cohorts). The selected model of variation in adult body mass included additive effects of absolute antler length (slope of 0.022 ± 0.08, N = 121) and population (Table

FIGURE 4

 4 FIGURE 4 Relationship between adult body mass (logtransformed) and fawn absolute antler length (log-transformed and measured at 8 months of age). Chizé: dark line, open circle; Trois-Fontaines: grey line, grey points. The relationship and data are represented on a double logarithmic scale

  t h e s e l e c te d m o d e l in cl u d e d b o d y m ass a n d a n inte r a c t io n b et we e n population and antler length (Table

FIGURE 5

 5 FIGURE 5 Relationship between the longevity and yearling relative antler length (at 20 months of age) for Chizé (dark line, open circle) and Trois-Fontaines (grey line, grey points). The relationship and data are represented on a double logarithmic scale

  though the low sample size at Trois-Fontaines prevents to draw any firm conclusion. Overall, yearling males that grow long antlers for their mass might play a fast life-history tactic, involving the repeated defence of large territories, but impaired survival. Short antlers in yearlings should delay territoriality establishment[START_REF] Vanpé | Age-specific variation in male breeding success of a territorial ungulate species, the European roe deer[END_REF], and thereby age at first mating. Delayed age at first

FIGURE 6

 6 FIGURE 6 Relationship between the probability of surviving beyond 6 years of age and yearling relative antler length (at 20 months of age) in Chizé (a, b) and Trois-Fontaines (c, d). Points (of different size according to the sample size) indicate the average survival probabilities within each class of yearling relative antler length, or body mass. In both cases, the interval between two successive classes corresponds to an increment of 5%

TABLE 1

 1 Summary of published studies that have evaluated the long-term life-history consequences of horn growth in male

bovids Species Population Horn growth Late-life trait Relationship References

	Bergeron, Festa-Bianchet, von Hardenberg,	and Bassano (2008)	Toïgo, Gaillard, and Loison (2013)		Geist (1966)			Bonenfant, Pelletier, Garel, and Bergeron	(2008)		Corlatti, Storch, Filli, and Anderwald (2017)						Loehr, Carey, Hoefs, Suhonen, and Ylonen	(2007)	Robinson, Pilkington, Clutton-Brock,	Pemberton, and Kruuk (2006)
	Null		Null (tend to be negative at very	old age)	Negative (greater horn growth in	the short-lived rams)		Null (once the negative bias due	to hunting mortality is removed)		Null		Negative		Null		Negative		Lower for N than S	No difference between N and S	Negative (N)	Null (S)	Null (N)	Null (S)
	Longevity		Age-specific survival		Survival at 7-11 vs.	12-17 (years of age)		Age-specific survival			Longevity						Longevity		Longevity	Lifetime breeding	success	Longevity	Lifetime breeding	success
	Early annual horn	increment	Second horn increment		Horn increment	between 3 and 8 years	of age	Horn increment	between 1 and 4 years	of age	Horn increment in the	first 2 years of horn	growth (L2 segment)				Horn volume from the	first 5 summers of	horn growth	Horn type: normal (N)	or scurred (S)	Horn size
	Belledonne-7 Laux Reserve,	France	Gran Paradiso National Park,	Italy	Banff National Park, Canada			Ram Mountain, Canada a			Swiss National Park,	Switzerland	Comprensorio Alpino di	Sondrio, Switzerland a	State Forest district of	Oberammergau, Germany a	Kluane National Park, Canada a		St Kilda, Scotland
	Alpine ibex (Capra ibex)				Bighorn sheep (Ovis	canadensis)					Chamois (Rupicapra	rupicapra)					Dall sheep (Ovis dalli)		Soay sheep (Ovis aries)

a Hunted populations.

TABLE 2 Sample

 2

		Independent variables
		Fawn antler	Yearling
	Dependent variables	length	antler length
	Longevity	137 (94; 43)	180 (101; 79)
	Survival beyond 6 years of	137 (94; 43)	190 (105; 85)
	age

sizes used per analysis. Samples sizes in brackets are population-specific (Chizé; Trois-Fontaines) Survival beyond 2 years of age 190 (105; 85) Adult body mass 121 (98; 23) 126 (78; 48) Adult antler length 111 (89; 22) 110 (70; 40) Antler length (20 months) 93 (76; 17)

TABLE 3

 3

	Parameter estimates from the model of variation in longevity as a
	function of relative allocation to antlers in
	yearlings (see Table S7b). Chizé is used as
	the reference population

TABLE 4

 4 Parameter estimates from the selected models of variation in the probability of surviving beyond six years of age (on a logit scale) as a function of relative allocation to antlers in yearlings (see TableS7b). Chizé is used as the reference population

		Estimates	SE	t
	Intercept	0.24	5.29	0.046
	Population	-11.06	5.20	-2.13*
	Standardized	-2.19	0.94	-2.33*
	yearling antler
	length
	Yearling body mass	3.46	1.32	2.61**
	Standardized	2.11	1.02	2.08*
	yearling antler
	length × Population
	*p < .05; **p < .01.

ACKNOWLEDGEMENTS

Field operations were supported by grants from the Agence Nationale de la Recherche (ANR-15-CE32-0002-01), and performed within the framework of the LABEX ECOFECT (ANR-11-LABX-0048) of Université de Lyon, within the program "Investissements d'Avenir" (ANR-11-IDEX-0007) operated by the French National Research Agency (ANR).

We thank all the ONCFS staff, in particular Claude Warnant and Gilles Capron, and the field volunteers for the organization of the roe deer captures. We warmly thank Dylan Childs and two anonymous referees for helpful and constructive comments that markedly improved our manuscript. The authors have declared that no conflict of interest exists.

The land manager of both sites, Office National des Forêts (ONF), permitted the study of the populations (Partnership Convention ONCFS-ONF, dated 23-12-2005). All procedures were approved by the Ethical Committee of Lyon 1 University (project DR2014-09, 5 June 2014).

Funding information ANR AGEX, Grant/Award Number: ANR-15-CE32-0002-01; LABEX ECOFECT, Grant/ Award Number: ANR-11-LABX-0048; ANR-11-IDEX-0007; French National Research Agency (ANR)

DATA ACCESSIBILITY

All data have been deposited on Dryad Digital Repository: https:// doi.org/10.5061/dryad.vt169tn [START_REF] Lemaître | Data from: The influence of early-life allocation to antlers on male performance during adulthood: Evidence from contrasted populations of a large herbivore[END_REF].

ORCID

Jean-François Lemaître

http://orcid.org/0000-0001-9898-2353

Louise Cheynel http://orcid.org/0000-0003-2292-4143

AUTHORS' CONTRIBUTIONS

J.F.L. and J.M.G. conceived the ideas and designed methodology; all authors collected the data on both study sites; and J.F.L., F.D. and J.M.G. analysed the data. J.F.L. wrote the first draft of the manuscript with input from J.M.G. and A.J.M.H., and all authors contributed critically to the drafts and gave final approval for publication.


Xiaohang Pan

Hao Xu

Yanqing Gao

Yafeng Zhang

Liaoxin Sun

Dan Li

Zhengji Wen

Shimin Li

Weiwei Yu

Zhiming Huang

Jianlu Wang

Bo Zhang

Yan Sun

Jinglan Sun

Xiangjian Meng

Xin Chen

Béatrice Dagens

Jiaming Hao
email: jiaming.hao@mail.sitp.ac.cn

Yue Shen

Ning Dai

Junhao Chu

Spatial and frequency selective plasmonic metasurface for long wavelength infrared spectral region

Keywords: Perfect absorber, long wavelength infrared, plasmonic metasurface, near-field microscopy, hyperspectral infrared imaging

The development of novel approaches that control absorption and emission operating in long wavelength infrared (LWIR) spectral region is of fundamental importance for many applications, such as, remote temperature sensing, environmental monitoring, thermal imaging, radiation cooling and industrial facility inspections. A high performance plasmonic metasurface-based absorber for the LWIR spectral region is presented. In our design, a pyroelectric thin film, poly(vinylidene fluoridetrifluoroethylene) (P(VDF-TrFE)) copolymer, is introduced as spacer, that offers the device not only with multiple selective high absorption bands but also promising potential for application in optoelectronics. By employing a scattering-type near-field optical microscopy (s-SNOM), both the near-field amplitude and phase optical responses of the absorber are investigated at resonant wavelength, thereby providing direct experimental evidence to verify the nature of the absorption effect. To further demonstrate the versatility of our design, a particular metasurface patterned by the building blocks of the plasmonic absorber is fabricated and characterized. Two-dimensional hyperspectral images show that such a patterned structure exhibits both frequency and spatially selective absorption.

Introduction

Long wavelength infrared (LWIR), in the commonly used division scheme of infrared detection, it corresponds to the electromagnetic radiation in the wavelength region roughly from 8 to 15 μm. [START_REF] Byrnes | Unexploded Ordnance Detection and Mitigation[END_REF] Most of the thermal radiation emitted by objects on earth at room temperature falls in the spectral range, then this region is also referred to as "thermal infrared". The development of novel approaches that control absorption and emission operating in this spectrum is fundamentally important for a variety of purposes, such as, remote temperature sensing, environmental monitoring, thermal imaging, radiation cooling, night vision and industrial facility inspections. [START_REF] Mahulikar | [END_REF][3][4][5] Metasurfaces are man-made ultrathin photonic materials constructed by subwavelength planar units with desired electromagnetic properties, which exhibit unusual abilities to manipulate electromagnetic waves, in particular, function as high-efficiency absorbing media. [6][7][8][9][10] In this letter, we propose and experimentally demonstrate a high performance plasmonic metasurface-based absorber (PMBA) operating in the LWIR spectral region. The basic structure of our design is comprised of a layer of metallic disk array separated from an opaque plasmonic ground plate by a dielectric spacer (See Figure 1(a)). It notes that such metal-insulator-metal (MIM) configuration metasurfaces have been suggested as absorbers and realized at various frequencies. [11][12][13][14][15][16] In contrast to previous designs, we incorporate a pyroelectric thin film, poly(vinylidene fluoridetrifluoroethylene) (P(VDF-TrFE)) copolymer (70/30 molar ratio), [17,18] instead of traditional dielectric materials (such as Al2O3, SiO2 and MgF2) for the dielectric spacer. This change not only brings our design with high performance multiple selective absorption bands resulting from the unique optical properties of the copolymer, but also makes it has great prospects for applications in compact optoelectronic devices such as uncooled LWIR pyroelectric detectors. [19][20][21] We note that perfect absorption in metasurface-based absorbers originated from localized electromagnetic resonances in the structures, the physical origin of the high absorption effect, has been pointed out. To illuminate the mechanism, electromagnetic field distributions at the working wavelengths obtained by numerical simulations have also been widely investigated so far. But there is little work to experimentally verify this argument. Here we map the spatial near-field distributions of our proposed absorber around resonance using phase-resolved scattering-type near-field optical microscopy (s-SNOM) and compare the experimental results with numerical simulations. [22,23] Finally, to further elucidate the functionality of our absorber design, we pattern a structure based on the PMBA, which exhibits both frequency and spatially selective absorption and investigate its corresponding hyperspectral thermal infrared imaging.

Results and Discussion

Because of the lack of the accurate optical constants of P(VDF-TrFE) copolymer, we firstly perform experiments to determine the optical properties of the copolymer by using infrared variable-angle spectroscopic ellipsometry (IR-VASE). Figures 1(c) and 1(d) respectively show the measured ellipsometry data Ψ and Δ (dashed-dotted curves) for a 200-nm-thick P(VDF-TrFE) copolymer thin film in the spectral range 600-6000 cm -1 (16.7-1.6μm) at incident angles θ = 60°, 70°. For definition of the ellipsometric parameters Ψ and Δ and detail information of spectroscopic ellipsometry, please refer to the ref. 24. The complex refractive index components n and k are then extracted by model fitting the experimental spectroscopic ellipsometric data (see Figures 1(c) and 1(d), solid curves), and presented Figures 1(e) and 1(f), respectively.

It is observed that the pyroelectric copolymer possesses a set of strong molecular vibrational modes in the LWIR spectral range. [START_REF] Bai | [END_REF][26][27] We now perform numerical simulations to investigate the absorbance (A) of the designed absorber with different values of the critical geometric parameters. The simulations are conducted using commercial package Lumerical FDTD solutions based on finite-difference-time-domain (FDTD) method. Owing to the presence of the continuous metallic ground plate (aluminum with the thickness h = 100 nm), the transmittance (T) of the device is nearly zero over the entire wavelength range we concerned. The absorbance (A) can thus be obtained using a simple formula A=1-R, where R denotes the reflectance, which can be directly obtained by numerical simulations. In the main text, we only concentrate on how the absorbance of the proposed device varies against the thickness (l) of spacer layers, the relationship between the absorption properties and the periodicity (p), the diameter (d) and the thickness (t) of metallic disks are presented in the supporting information as shown in Figures S1, S2 and S3, respectively. [28] Figure 2(a) shows the simulated absorbance spectra for four different thickness spacer layers (l = 0.3, 0.37, 0.45, and 0.52μm) at normal incidence, with other parameters d = 3 m, t = 50 nm, and p = 6 m. As shown in Figure 2(a), when thickness l = 0.37μm (denoted by "PMBA-S1" hereinafter), the absorbance can be as high as 99.5% at the wavenumber of 980 cm -1 (10.2μm). For l = 0.52μm case, at the wavenumber of 825 cm -1 (12.1μm), the maximum absorbance reaches over 99% as well.

In order to experimentally demonstrate this effect, the absorbers are fabricated on a polished silicon wafer with standard microfabrication techniques. Figure 1(b) displays a scanning electron microscopy (SEM) image of a fabricated sample which adopts the optimized geometric parameters of PMBA-S1. Transmission and reflection of these devices are characterized using a Fourier transform infrared (FTIR) spectrometer (Bruker IFS 125HR). The transmission is measured at the normal incidence and reflection is measured at incident angle of 15°. Figure 2(b) shows the spectra of experimental absorbance, which are obtained based on the measured transmittance and reflectance. Reasonable agreements are found when compare the experimental spectra with the numerical simulated. For thickness l = 0.37 (0.52) μm, the maximum absorption is close to 95.5% (97%) at the wavenumber of 981 (829) cm -1 in experiment. In addition, from both the numerical simulated and experimental results, one can see that for each case, the absorber exhibits three distinct absorption peaks. More notably, it is found that all of these high performance absorption peaks are attributed to the excitation of the fundamental cavity resonant mode (The evidence is presented in the supporting information as shown in Figure S4). This is in sharp contrast to the multispectral metamaterial absorbers reported in the literatures, in which the multiband high performance absorption is due to either the hybrid multi-resonant unit structures or simultaneously exciting several high-order localized electromagnetic resonant modes. [29][30][31][32] We further reveal that such multiband absorption effect is strongly related to the unique optical constants of the P(VDF-TrFE) copolymer spacer. To figure out the physical phenomenon, time-averaged dissipative power (TADP) is theoretically investigated. For example, Figure 2(c) shows the normalized TADPs (Ploss / Pin) spectra for the total and each component of the PMBA-S1. As can be observed, for the two absorption peaks of 825 cm -1 (12.1 μm) and 869 cm -1 (11.5 μm), more than 78% and 87% energy of incident radiation are dissipated by the copolymer spacer, respectively, while for the absorption peak of 980 cm -1 , energy dissipation by the spacer is less than 58%. This is exactly thanks to the optical properties of the copolymer spacer, which exhibit two strong infrared-active vibrational modes around the two long-wavelength absorption peaks and have rather larger attenuation factors compared to that of the short-wavelength absorption peak (see Figure 1(e)).

To visualize and better understand the nature of our absorber device, we utilize s-SNOM (Neaspec) to investigate the near-field response of our fabricated absorber, which can provide near-field optical signal amplitude and phase by interferometric detection of scattered light with a phase-modulated reference wave. Figures 3(a Numerical simulations are also conducted for this structure. The corresponding amplitude and phase near-field images are displayed in Figures 3(c) and 3(d). Good agreements are noted between the numerical and experimental results.

In order to investigate the angular response of optical performance of the PMBAs, additional optical measurements and numerical simulations are performed. show that for TM polarization, a clear anticrossing behavior is observed for the absorption peak 980 cm -1 at the incident angle of 45    . The reason for this angular characteristic is that the oblique incident wave is diffracted by the first-order reciprocal vectors of the periodic disk array, and exciting the surface plasmon polariton (SPP) supported by the air-spacer-metal multilayer, which thus strongly couples with and disturbs the localized resonant mode. [33][34][35] Nevertheless, for the TE polarization, the anticrossing point is absent. This is due to that although the diffracted wave is also introduced for high incident angles, which is no longer interfering with the localized resonant mode, since the multilayer cannot support TE surface modes. [START_REF] Maier | Plasmonics: Fundamentals and Applications[END_REF] It is worth mentioning that diffraction beams are introduced by the periodic structures for both polarizations and play the roles of new channels for wave scattering at high incident angles. As a result, the absorption performances would be highly decreased.

This opens a window for our design to function as optical devices with directional and frequency-selective absorption/emission properties. [35] For further demonstration of the versatility of the PMBA and its prospect in practical applications, we experimentally pattern a metasurface, which exhibits both frequency and spatially selective absorption properties. [START_REF] Liu | [END_REF] Figure 5(a) shows an optical microscope image of the patterned metasurface, printing the two letters "IR"

(the abbreviation of Infrared). In this design, the basic building block of the plasmonic absorber PMBA-S1 is used to form the pattern of the letters. The pattern contains a total of more than 22000 plasmonic unit cells. The minimum feature size of the letters is about 150μm. The unit cell of background area has the identical lateral dimensions as S1, i.e., 6μm × 6μm, and consists of the same metallic ground plate and P(VDF-TrFE) copolymer with identical thicknesses, but no metallic disk on the top.

We characterize the sample using a FTIR combined with an infrared microscope. The sample stage can be moved automatically and is controlled by a computer. A spatial resolution of 50μm is adopted, thus generating an image with 40 pixels along both x and y axes. For each pixel, the reflection spectrum is collected in the ranges 600 -6000 cm -1 with the spectral resolution of 2 cm -1 . Eventually, two-dimensional hyperspectral images are obtained, with a total data record size of 4.3×10 6 . These allow us to fully study the absorption property of the sample.

The experimental absorbance spectrum of the background is shown in Figure 5

Conclusion

In summary, we have presented a plasmonic metasurface-based ideal absorber for the LWIR spectral range. In our design, a pyroelectric thin film, P(VDF-TrFE) copolymer, is introduced as the spacer layer, which not only brings our device with multiple selective high absorption bands, but also offers a further step towards the development of the compact optoelectronic devices. The angle resolved optical responses show that the absorber is sensitive to the incident angles and can be controlled by the periodicity of the disk array. This opens another window for our design to function as optical devices with directional and frequency-selective absorption/emission characteristics. To better understand the origin of the absorption effect, we applied interferometric near-field microscopy to map both the near-field amplitude and phase optical responses of the absorber, and the obtained images provide direct experimental evidence to confirm that such high performance absorption effect is exactly due to the excitation of the fundamental electromagnetic resonant mode. We further patterned a metasurface based on the PMBA, which exhibits both frequency and spatially selective absorption. The above listed advances suggest that the absorber has various promising potential applications, including thermal infrared light sources, thermal management, and hyperspectral thermal imaging. Optical measurements. The angle resolved reflection spectra measurements are performed by using a Fourier transform infrared spectrometer (Bruker IFS 125HR), which is equipped with reflection module allowing for angles ranging for 15° to 80°.

The measured reflection spectra are normalized with respect to a gold mirror. A wire-grid infrared polarizer is placed in the light path to obtain both TE and TM polarization spectra. Two-dimensional hyperspectral images of the pattern sample are characterized with the Thermo Fisher microscopic Fourier infrared spectrometer.

Near field microscopy. The near-field optical response amplitude and phase of the absorber devices are characterized by employing a commercial s-SNOM (Neaspec GmbH, Germany). Our near field microscopy is equipped with a tunable quantum cascade laser with operating wavelengths range of 7.8 -10.5 m. Nano-scale-resolved images can be obtained by collecting of the infrared light scattered at an atomic force microscope tip with a phase-modulated reference wave.

) and 3 (

 3 b) show the experimental near-field amplitude and phase images of vertical electric field component Ez for the l = 0.37 μm structure recorded at λ = 10.2 μm. As one can see, the amplitude shows noticeable enhancement on the surface of metallic disks, the phase signal has a remarkable variation of about 180° at the disk center. The images provide direct experimental evidence that such high performance absorption effect is as expected owing to the fundamental dipolar-like near-field resonant mode.

Figures 4 (

 4 Figures 4(a) and 4(b) display the measured reflectance of the sample shown in Figure

  (b) as the dashed blue curve (denoted by Area 1), which is recorded at the position labeled by the blue spot in Figure5(a) and exhibits low absorption performance over the LWIR spectral range. For comparison, in Figure5(b) we also plot the absorbance spectra of the area composed of PMBA-S1 (solid red curve, Area 2) recorded at the position labeled by the red spot in Figure5(a). As expected, the absorption of this area is identical to the one presented in Figure2(b) as the solid red curve. The 2D and 3D images of the experimental absorbance data for the absorption peak wavelength of 10.2 μm (980 cm -1 ) are displayed in Figures5(c) and 5(d), respectively. As can be observed from these two images, at this resonant wavelength, the pattern composed of the PMBA achieves nearly ideal absorption. To illustrate the frequency dependence of the spatially selective structure, Figures5(e) and 5(f) respectively show 2D and 3 D images for a wavelength of 14 μm (714 cm -1 ). We can see that in contrast to the optical photograph (Figure5(a)) and the images at 10.2 μm (Figures5(c) and 5(d)), at 14 μm the pattern cannot be so clearly distinguished from the background.

4 .

 4 Experimental section: IR-VASE characterizations. A 200-nm-thick P(VDF-TrFE) copolymer (70/30 molar ratio) thin films are prepared on electronic-grade Si(100) wafer by spin-coating method for ellipsometric measurements. Spectroscopic ellipsometry measurements are performed using an infrared variable angle spectroscopic ellipsometry system (IR-VASE, J. A. Woollam). Lorentz model with multiple oscillators is used to fit the experimental data. Numerical simulations. Full field electromagnetic wave simulations are conducted using commercial package Lumerical FDTD solutions based on finite-difference-time-domain method. In our simulations, both the optical constants of gold and aluminum are taken from Ref. 38. Reflection and Transmission are respectively collected with power monitors behind the incident plane wave and the structure. Electromagnetic field distributions are recorded by 2D field profile monitors at the corresponding observation positions presented in the text. Sample fabrication: The LWIR absorber samples are fabricated using the techniques of direct laser lithography and wet-chemical lift-off. Clean silicon substrates are coated with the bottom Al layer and P(VDF-TrFE) copolymer spacer by magnetron sputtering and spin-coating methods, respectively. And then, the disk array is defined by in AZ5214 photoresist using direct laser lithography. A 50-nm-thick gold is then deposited using magnetron sputtering. Complete devices are achieved by removing the photoresist mark in alcohol solutions.

Figure 1 .

 1 Figure 1. (a) Schematics of the LWIR metasurface-based absorber. (b) SEM image of an experimental sample. (c) and (d) IR-VASE data of Ψ and Δ (dashed-dotted curves) for a 200-nm-thick P(VDF-TrFE) copolymer thin film at angles of incidence θ = 60° and 70°, and the fitted results from the multiple oscillators Lorentz model (solid curves). (e) Real and (f) imaginary parts of refractive index for the P(VDF-TrFE) copolymer obtained from the Ellipsometry data analysis.

Figure 2 .

 2 Figure 2. (a) Simulated and (b) measured absorbance spectra for the samples with four different thickness P(VDF-TrFE) copolymer spacers. (c) The normalized time-averaged dissipative power spectra for the total and each component (top Au disk, P(VDF-TrFE) copolymer spacer, and bottom Al mirror) of the l = 0.37 μm absorber device.

Figure 3 .

 3 Figure 3. (a) and (b) Experimental near-field amplitude and phase images of the vertical electric component for the l = 0.37 μm absorber device mapped at the wavelength of 10.2 μm. (c) and (d) Simulated near-field amplitude and phase images of the corresponding model system.

Figure 4 .

 4 Figure 4. (a) and (b) Measured reflectance of the l = 0.37 μm absorber device as a function of wavenumber and in plane wave vector kx for both TE and TM polarizations. (c) and (d) The corresponding simulated reflectance spectra for both TE and TM polarizations.

Figure 5 .

 5 Figure 5. (a) An optical microscope image of the patterned spatially dispersive metasurface. (b) The experimental absorbance of the background and the area composed of the metasurface-based absorber. (c) 2D and (d) 3D images of the experimental absorption at 10.2 μm. (e) 2D and (f) 3 images of absorption at 14 μm.

Acknowledgements

This work is supported by National Key R&D Program of China (2017YFA0205801), National Natural Science Foundation of China (61471345), the National Young 1000 Talent Plan, GDRI Sino-French Optoelectronics and Photonics -Photonet, and Shanghai Science and Technology Committee (16JC1403500).