Willy Cottin
email: willy.cottin@stellantis.com

Yuqi Liu

Guillaume Colin
email: guillaume.colin@univ-orleans.fr

Alain Charlet

Sébastien Houillé

From boundary lines to boundary surfaces for dynamic programming with final state constraints ⋆

Keywords: Dynamic Programming (DP), Boundary Lines (BL), Boundary Surface (BdS), Optimal Control Problem (OCP), Energy Management Strategy (EMS)

à la diffusion de documents scientifiques 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.

Introduction

Vehicle powertrain electrification is considered an efficient solution to the problem of reducing CO 2 emissions from the automotive industry. Electrification requires multiple power sources in one vehicle, making a control to distribute the powers, called Energy Management Strategy (EMS), necessary. To evaluate the performance of strategies, criteria are mandatory. Moreover, an optimal control policy can help to correctly size a given system. Naturally in the automotive field, the first criterion to be minimized is the fuel consumption. However, with electrification, it is more appropriate to minimize energy consumption. These strategies can be written as Optimal Control Problems (OCP) and solved using optimal control methods. One method is Dynamic Programming (DP), a numerical method developed by R. Bellman in the 50s [START_REF] Bellman | The theory of dynamic programming[END_REF]. The OCP can also be written in its dual form and solved with the Pontryagin Maximum Principle (PMP) [START_REF] Pontryagin | Mathematical theory of optimal processes[END_REF].

Since PMP is an analytical method, the problem resolution is faster than with DP resolution. However, the main drawbacks of this method are constraints and non linearity considerations throughout the simulation, which cannot be easily handled with PMP but are straightforward with DP. This is mainly why DP is one of the reference methods to solve OCP. As DP is a numerical method, it explores all the possibilities and then determines the optimal control ⋆ This paper is supported by the Association Nationale de la Recherche et de la Technologie (ANRT) through a CIFRE contract. It is also developed under the collaborative framework OpenLab Energetics between Stellantis and PRISME laboratory.

policy to apply to the dynamic system while minimizing a given criterion.

Because of the computation time and memory complexity of DP, research has been conducted to optimize it in order to handle increasingly complex systems. For example, Approximate Dynamic Programming (ADP) [START_REF] Bertsekas | Approximate dynamic programming[END_REF] replaces the cost-to-go matrix in the backward step to avoid multiple interpolations. A second possibility is to use Iterative Dynamic Programming (IDP) [START_REF] Luus | Iterative dynamic programming[END_REF]. This solution aims to reduce time calculation by using a coarse mesh grid for all domains. Then, the same mesh grid is reused but domains are contracted around the previous solution. This provides an interesting accuracy for a reasonable calculation time. Another branch of research to reduce calculation time and calculation complexity has focused on Boundary Lines (BL) [START_REF] Sundström | On implementation of dynamic programming for optimal control problems with final state constraints[END_REF]. This consists in determining the extreme states that can be visited while satisfying final constraints. In the literature, a single state boundary line problem is well known. However the resolution of boundary lines for two states or more, called boundary surface (BdS) here, does not have an explicit solution. The study in [START_REF] Van Schijndel | Dynamic programming for integrated emission management in diesel engines[END_REF] attempted to estimate BdS for a two-states problem in order to mesh only in this area. This produces a more accurate result with the same mesh resolution.

In section 2, the Optimal Control Problem formulation is provided. Then in section 3, Dynamic Programming and the boundary surface method are explained. These methods are applied on two models in section 4. Section 5 concludes the paper.

The Optimal Control Problem formulation for a system consists in minimizing a given criterion while respecting constraints. The cost function J(x, u) is composed of the instantaneous cost L (x(t), u(t)) that is the consequence of the inputs applied to the system and states which the system passes through. The second element Φ (x(t f )) is a penalty term applied as a function of the final system state. The OCP in its continuous form is the following one:

           min u J(x, u) = t f t 0 L(x(t), u(t))dt + Φ x t f w.r.t.: ẋ(t) = f (x(t), u(t)) x(t) ∈ Xt, ∀t ≥ t 0 u(t) ∈ U t, ∀t ≥ t 0 (1)
where x and u are respectively the state and the input vectors. ẋ is the state derivative (the dynamics) which depends on the inputs and the current states. The input and state sets are defined as:

U = {u ∈ R n , u min (t) ≤ u(t) ≤ u max (t)} X = {x ∈ R m , x min (t) ≤ x(t) ≤ x max (t)} (2)
Since DP is a numerical method, the OCP will be written in discrete form which requires a grid of inputs, states and time, giving:

           min u J(x, u) = N -1 k=0 L(k, x k , u k ) + Φ (x N ) w.r.t.: x k+1 = f k (x k , u k ) x k ∈ X k , ∀k ≥ 0 u k ∈ U k , ∀k ≥ 0 (3)
Here, Φ (x N ) is evaluated at each state mesh point at final time. There are several ways to treat this. One is to impose an infinite penalty for all final states that do not respect the constraint, as shown in equation (4).

if

x N ∈ X N , Φ(x N ) = 0 else, Φ(x N ) = Inf (4)
Instead of considering infinity for Φ(x N ) , it is also possible to have a sufficiently large number which will help the algorithm to find solutions. Another possibility is to consider some smooth penalty. This is done by using a given penalty function depending on final states. This approach will not lead to an infeasible solution but it has to be well tuned to converge to the target states.

3 Dynamic Programming

Dynamic Programming algorithm

Dynamic Programming is a numerical method to solve OCPs. To do so, the discrete dynamic programming principle is used:

V (k, x k ) = min u (L(k, x k , u k ) + V (k + 1, x k+1 )) (5)
where V represents the cost-to-go function evaluated at each point in the space domain. The algorithm is divided in two parts. The first one, called backward, is a time loop going from the end to the beginning of the simulation. During this loop, the model is run in forward form:

x k+1 = f k (x k , u k ) with
x k all the possible states and u k the possible inputs according to the mesh grid. Then the instantaneous cost L(k, x k , u k ) is determined. However, to compute the whole cost-to-go matrix, an interpolation is mandatory to estimate the value for all possibilities. This interpolation can be either a nearest one, linear or any other interpolation method. Finally, equation ( 5) is solved for every state mesh node and the inputs associated to this optimal cost are stored in an optimal input matrix.

The second part, called forward, is also a time loop but going from the beginning to the end of the problem, where optimal inputs are determined by using results from the previous part. Moreover, the optimal policy is applied and the optimal cost is calculated at each time step. In the end, the whole system trajectory is determined as well as the global optimal cost.

To respect the final constraints, the optimal cost-to-go matrix is initialized with some values. These values were explained in the previous section. If a hard penalty solution is used, some interpolation errors may occur due to the large values of penalties. Knowledge of the boundaries reduce/avoid contamination of the cost-to-go and improve DP accuracy.

Boundary Lines (BL)

Because the cost-to-go matrix can be contaminated with penalties, this may lead to a poor construction of the costto-go matrix and furthermore to a poor result. The role of boundary lines is to determine exactly where the frontier between the feasible and infeasible domain is. Thus, BL prevents this contamination from occurring. BL operate with some assumptions. A first assumption is that all possible states existing in the boundary area are really feasible. This means that if a state is included in the area described by the boundary contour, then there is necessarily a solution to achieve the simulation in the final boundary area. In other words, it is impossible to have an area of infeasible points inside the boundary area. A second assumption is that the model used can be inverted.

If not, it is impossible to do step 2.(a) of the following algorithm. BL are then determined by using the dynamic equations in the inverted form

x k = x k+1 -f -1 k (x k+1 , u k ).
The algorithm used to determine BL is as follows: 1.Initialize the x end with target boundary final states. 2.Backward loop in time for k = N -1 to 0: (a) determine upper and lower boundary lines at step k from the inverted dynamic equations and the upper/lower positions at time k + 1. (b) store state solutions, associated input(s) and cost for upper and lower boundary. If the state solutions at time k are outside the constraints, it is possible to enforce boundaries on these constraints and determine the associated inputs and costs. Lastly to determine boundaries and use the results during the backward step of DP, another time loop can be done before the backward one or it can be directly merged in the backward loop.

An example of the use of boundaries is shown in blue in Fig. 1. From the point x i+1 k+1 it is possible to reach the blue square which has no cost-to-go information. In order to estimate this value, an interpolation with the two nearest meshed points is done. In this example, one of the two points is the lower boundary line. This interpolation can be linear, nearest, etc. Thus the cost-to-go at the blue square with a linear interpolation is:

V (k + 2, x blue k+2 ) = V (k + 2, x i+2 k+2 ) + (x blue k+2 -x i+2 k+2 ). V (k, x bl l k+2 ) -V (k + 2, x i+2 k+2 ) x bl l k+2 -x i+2 k+2 (6) X i+3 J k i+3 J k+1 i+3 J k+2 i+3 X i+2 X i+1 X i k k+1 k+2 Time index k Infeasible Feasible State variable x J k i+2 J k+1 i+2 J k+2 i+2 J k i+1 J k low J k+1 low J k+2 low J k+1 i+1 Inf Inf Inf Inf
Fig. 1. Dynamic programming with cost-to-go of mesh and boundary lines scheme Fig. 1 shows the cost-to-go of meshed points (in space and time) whether in the feasible domain or not. Here, points outside the feasible zone have a cost-to-go set to infinite. Those zones are determined with boundary lines (the lower one in red here). In this scheme, all meshed points below the lower state at a certain time index are considered to be infeasible. Moreover, for each point belonging to the boundary line the cost-to-go is known and has been used to determine the cost-to-go for meshed points in the feasible domain. Once more, DP without BL cannot find a solution at the border of the interval because of the proximity to infeasible zones.

Boundary Surfaces (BdS)

To generalize BL to more complex systems, borders have to be determined whatever the problem dimension is. However, when the system state dimensions are higher than one, determining the boundary borders is more difficult than for a single state problem. Because of the complexity of the system, the upper/lower boundary situation cannot be repeated since the borders are now continuous for a time step and not two points, as is the case with boundary lines.

In [START_REF] Van Schijndel | Dynamic programming for integrated emission management in diesel engines[END_REF], the reachable set of states is determined from the union of the previous backward reachable set and the convex set defined by the maximum estimated values of each state. This method gives an approximation of the reachable set of states which is used to refocus the meshing only in this space. Even if some infeasible spots exist in the reachable set the mesh resolution will naturally exclude them because of their results. However, because of the curse of dimensionality of DP, a solution that works with a coarse mesh is interesting. Thus in this article, the idea is to determine the reachable set of states using only the previous meshed set in backward. This makes it to work only with the initial mesh (coarse of fine) without any change over time.

Continuity of the borders involves the estimation of all possibilities from the boundary surface to determine a new boundary surface. The boundary surface method is used numerically such that non-analytical problems can be solved. The boundary X k is the set of points at the edge of the feasible and infeasible zones. First of all, from the previous boundary set, all the inputs are applied to find all the points leading to the boundary P k .

P k = f -1 (X k+1 , u k ) (7)
In this set, some points are outside the constraints O k so they are removed Q k = P k \O k . Then, only points that are on the edge are kept to build the new boundary set (8). Qk is the closure while

• Q k is the interior of the Q k space. X k = Qk \ • Q k (8)
This is done from k = N -1 to 1 as in backward dynamic programming.

To limit the calculation time and memory resources required for boundary surfaces, it is possible to sample points belonging to this boundary. However, this sampling has to be done to conserve as much information as possible since points on the boundary could be massively placed in specific areas while some other areas could have very few points of information. Thus, if areas with fewer points are removed by the sampling, then there is a huge loss of information whereas if many points are removed from the massively populated areas, the loss of information is negligible. To perform the sampling, the perimeter of the boundary is divided into multiple parts. If the number of points in one part is higher than the maximum number of points for a given part, then a sampling is done in between all the points existing in this precise part. With this solution, only parts with an excessive number of points are reduced to the maximum number of points available and parts with few points will keep all their points.

The boundary surface was determined here by using the boundary Matlab function. This function takes as input coordinates of the studied points along each dimension and returns indexes of points belonging to the boundary surface.

Fig. 3 shows an example of an optimal trajectory with boundary surfaces on the same figure. The first model used to illustrate the benefits of using BdS is a simple convex problem with one command and two states. This model ( 9) considers a vehicle trying to reach a fixed speed v f and position p f at a precise timing t f . Thus the system input is the acceleration u(t) = a(t) and states are speed v(t) and position p(t).

ẋ = v ṗ = u v (9)
Here, the vehicle has to reach the given position p f = 200 m and speed v f = 100 km/h at a given time t f = 10 s. The optimization criteria is the integral of the square of the instantaneous acceleration, which leads to : L(x(t), u(t)) = a(t) 2 .

Optimal Control Problem resolution with Pontryagin Maximum Principle

Since the problem is convex, it is possible to use PMP to solve the OCP and have a reference to compare DP with several solutions: with or without BdS and with smooth penalty.

In this case the Hamiltonian H can be written as:

H(x(t), u(t), t) = L(x(t), u(t), t) + λ T ẋ =⇒ H(x(t), u(t), t) = u 2 (t) + λv(t).u(t) + λp(t).v(t) (10)
where λ p (t) and λ v (t) are the co-states to be found.

As H is convex, the optimal input u * (t) can be found as follows:

∂H(x(t), u(t), t) ∂u * (t) = 0 =⇒ u * (t) = -λ v (t) 2 (11)
while the co-state dynamics are defined by Eq. ( 12). λi (t) = -∂(H(x(t), u(t), t)) ∂x i (12) w.r.t.:

t f t0 u * (t)dt = v f -v 0 , t f t0 v * (t)dt = p f -p 0 (13)
From ( 13), it is possible to have analytical expressions of the two co-states λ v (t) and λ p (t) and then include them in expressions of optimal acceleration, speed and position to fully resolve the problem.

Optimal Control Problem resolution with Dynamic Programming

Dynamic programming is also used to solve this simplified acceleration problem with constraints on speed and position. The objective is to compare DP with and without the BdS and also DP with smooth penalty (DPSP) to PMP. ,2) = 500. This smooth penalty aims to have both states as close as possible to their targets. To achieve this coefficients were tuned by hand. Results show that DP without BdS has the highest cost and is the farthest away from the final constraints. DPSP satisfies the constraints and has a much lower cost than DP without BdS. However, the best DP solution is the one with BdS.

(x N ) = 10 20 if x N / ∈ X N otherwise, Φ (x N ) = 0. For the smooth penalty DP, Φ (x N ) = 2 i=1 C i .(x i N -x i f ) with x 1 = v, x 2 = p and C (1
All solutions are shown in Fig. 4 where the red dashed curve is the PMP solution, the blue dashed curve is DP without BdS, the green dashed curve is DP with smooth penalty and the black curve is DP with BdS. All accelerations have the same global behavior, i.e. full acceleration at the beginning followed by a linear decrease until the end. The PMP gives a pure linear acceleration while the three DP solutions show some fluctuations on the input linked to the numerical resolution of the problem. It is also important to note that the fluctuations are higher and occur throughout the simulation for DP without BdS, which leads to a sub optimal solution. The behavior of DPSP is close to that of DP with BdS. However, fluctuations at the end are higher, leading also to a sub-optimality compared to DP with BdS. Thanks to boundary surfaces, fluctuations are lower and appear at the end of the simulation. This behavior appears because the solution lies on or close to the boundaries. Nonetheless the fluctuations remains in the feasible domain throughout the simulation. At the beginning in time of the boundary surface DP, the acceleration is much smoother because of the use of boundary surfaces during the backward part to improve DP accuracy. The second model used to apply the methodology is a non convex, non linear model which corresponds to a parallel hybrid electric vehicle [START_REF] Jager | Optimal control of hybrid vehicles[END_REF]. Here, the vehicle powertrain is composed of a battery, a motor and an irreversible power source (engine, turbine or other). The irreversible power source is assumed to give only the exact required power instantaneously. To calculate the required power from the powertrain to follow the cycle, Newton's second law is used. Fig. 5 shows a block scheme of the studied system. To propel the vehicle at a certain speed, the powertrain has to provide wheel power P W . This power can be provided by the engine P E and/or the motor P M . To drive the motor power, the torque motor C M is used.

The battery State Of Charge (SOC B ) is the ratio of the actual capacity compared to the maximum capacity of the battery. Moreover, to regulate the motor temperature T M with the coolant mass flow ṁc , the battery provides power P P to the coolant system. The battery also provides power to the motor to ensure the motor power. The battery has to ensure the sum of these two powers P K expressed in eq.( 14).

       P K = η M .P M if P M < 0 P K = P M η M else with: η M = C M (t).ω M (t) C M (t).ω M (t)+R• C M (t) k 2 (14)
For more information on the modeling, see [START_REF] Cottin | On the use of frequency analysis tools to minimize the calculation complexity of optimal control problems[END_REF]. The 2 inputs 2 outputs system can be modeled as follows:

x = SOC B T M u = C M ṁc (
ṠOC B = (θ • ṁc(t) • Cpc • (T M (t) -Tc) + R • C M (t) k 2 + C M (t) • ω M (t)). -1 Q • (E + e • SOC B (t)) (16) ṪM = R • C M (t) k 2 -ṁc (t) • Cp c • (T M (t) -T c ) m m • Cp m (17)
Eq. ( 16) describes the battery State Of Charge SOC B variation over time. It depends on the power spent to use the cooling system θ • ṁc (t)

• Cp c • (T M (t) -T c ), the motor losses R • ( C M (t) k ) 2 and the propulsion C M (t) • ω M (t)
. The first power varies over time since it depends on the mass flow input ṁc (t) and the motor temperature state T M (t). The cooling temperature T c , thermal capacity Cp c and the thermal to electric power coefficient θ are constants. Power losses are defined by the constant motor resistance R and its intensity which is function of the torque motor input C M (t). The required power to propel the vehicle is a function of the motor torque and speedω M (t). These two variables vary in time. In order to have an intensity, electrical powers are divided by the battery voltage which is modeled here by a constant voltage E plus a voltage function of the battery State Of Charge e • SOC B (t). Thus this voltage is a function of time. Lastly, this intensity is divided by the maximum capacity of the battery Q. Eq. ( 17) which captures the motor temperature dynamic over time, is a thermal power balance. It describes the dynamic function of the heating power, losses here, and the power removed by the cooling system.

In this optimal control problem, the objective is to minimize the overall vehicle energy consumption while managing the battery state of charge and the motor temperature. The trip to achieve here is the WLTC cycle. The variables' constraints are expressed in (18).

     0 ≤ SOC B (t) ≤ 1 20 ≤ T M (t) ≤ 60 -100 ≤ C M (t) ≤ 100 0 ≤ ṁc(t) ≤ 0.05 (18)
The instantaneous criterion is expressed in (19). The objective here is to ensure a trip repeatability while optimizing the motor temperature.

L(x(t), u(t)) = P 2 E (u(t)) + f (T M ) • P 2 B (u(t)) with f (T M ) = aT 2 M + bT M + c (19)
Coefficients a, b and c are representative of the electrical efficiency of the powertrain. Fig. 6. Parallel hybrid electric vehicle optimal solution with boundary surface Fig. 6 presents the solutions given by DP without BdS with 100 mesh points for each variable in dashed red and with mesh points for variable in dashed black.

BdS is also shown with blue lines while the BdS projected in one state space and time is presented as a dashed blue line. The solution found by DP with BdS uses the mass flow coolant during all the trip compared to DPs without BdS. Thus the motor temperature is maintained at the right temperature with respect to the constraints. Moreover, the distribution of powertrain power seems to be more efficient to achieve the trip at the minimum acceptable final SOC B . The results is a lower energy consumption since there is no unnecessary battery filling.

Even with the small number of meshes for the BdS, the management of the final constraints SOC B ≥ 0.61 and T M ≥ 37°C is much better. It can be seen that BdS is closer than classical DP to the minimum final constraint value, leading to a reduced sub-optimality. The results are shown in Table 2 The above results show the final system states for three different simulations shown in Fig. 6. The computer used is an Intel CORE i7 processor at 3.6 GHz with 96 Go of RAM. It shows the utility of the BdS through the expected cost, which is the smallest one despite the allocation of resources. However, for the same mesh grid, the total calculation time is much higher because of the interpolation function while considering the boundary surface points. To achieve the same optimal result, BdS required less time than DP without BdS, making it a valuable option. On top of this, it might be possible to reduce this calculation time and thus obtain even greater benefit from the BdS.

Conclusion

This article presents a solution to extend the Boundary Lines method for one-state problems to a Boundary Surface method for two-states problems. Usage has been shown on two problems: one convex problem solved with PMP and DP and one non-convex problem solved only with DP with/without the BdS. The idea was to estimate the usefulness of the BdS for multi-state systems in the resolution of the optimal control problem.

For a simplified convex problem, a comparison between PMP and DP with/without BdS and DP with smooth penalty was done. It shows that for the same mesh grid, if the smooth penalty is well tuned, BDSP can satisfy final constraints with a lower cost than DP without BdS. However, BdS offers results much closer to the reference set by PMP than DP without BdS or with smooth penalty.

The boundary surface solution was then tested on a convex and a non-linear and non convex system to compare the results of DP with and without BdS. In this case, a comparison with PMP is no longer possible. Results show that the calculation time is much higher with BdS for the same mesh grid but that the gain in accuracy is very high. In such that it is possible to have a better solution than the finest mesh grid possible with an average computer for a reasonable calculation time.

Based on the way this BdS works, it is possible to use the proposed method on systems with more than 2 states. For systems with 3 states, the boundary will be edges of a volume and for more than 3 states, the feasible area will be an hypervolume. Since the calculation time with BdS can be high, future work will focus on reducing the time needed for the BdS method.

Fig. 3 .

 3 Fig. 3. BdS illustration with optimal trajectory in state space over time

Fig. 4 .

 4 Fig. 4. Acceleration problem resolution with PMP, DP without BdS and DP with BdS

Fig. 5 .

 5 Fig. 5. Parallel Hybrid Electric Vehicle block scheme with the following dynamics:

  DP without BL arrives in the center of the target interval while DP with BL is led by the lower BL.

		1			Interest of Boundary Lines
				DP with BL trajectory
		0.9		DP without BL trajectory
				Boundary lines
		0.8
		0.7
	State (-)	0.5 0.6
		0.4
		0.3
		0.2
		0.1
		0
		0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1
							Time (-)
	Fig. 2. Illustration of the interest of the boundary lines
	method
	Fig. 2 shows the trajectory difference between a classic
	DP in a blue dashed line and DP with BL in black.Here
	the boundary lines are shown in red dashed lines and the
	same mesh is used for both DPs. At the beginning, the

trajectories are similar. However, from the first half of the simulation,the black line drops below the blue one. The division appears because of the boundary method. Thanks to it, the cost-to-go matrix is now healthy in the area close to the lower boundary lines. This leads to a better result by decreasing the state in this example. At the end of the simulation,

Table 1 .

 1

				DP w/o	DP w/
	Methods	PMP	DPSP	BdS	BdS
	Accel. Mesh (-)	-	201	201	201
	Speed Mesh (-)	-	201	201	201
	Position Mesh (-)	-	201	201	201
	Time step (s)	-	0.1	0.1	0.1
	Cost (-)	122	128.5	175.5	124.5
	Final Speed (km/h)	100	100	110	101
	Final position(m)	200	200	208	200

Table 1 sums up the main results. In Table 1, penalty Comparison of results between PMP and DP with/without BdS and DPSP Φ

  .

		DP w/o	DP w/o	DP w/
	Methods	BdS	BdS	BdS
	Mesh points (-)	40	140	40
	Cost (MJ)	11.1	8.3	6.8
	T M (t f ) (°C)	40.2	37.5	37.4
	SOC B (t f ) (-)	0.78	0.7	0.61
	Calc. Time (s)	336	55 168	55 566
	Table 2. Result comparison between DP with-
		out and with BdS


François Warlop
email: francois.warlop@grab.fr

Melodie Aujogue

Alice Bombeau

Mathias Boucheraki

Lou Milon

Raphael Paut

A user-friendly decision tool for facilitation of fruit trees choice

Keywords: decision tool, fruit trees, design decision tool, fruit trees, design

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

In France, farmers are getting more and more interested in integrating trees in their agroecosystems, thanks to the French plan called « Plantons des haies » giving substantial subsidies for plantations. Breeders also consider the multiple interests of trees for broilers or even ruminants. Fruit trees are among the preferred tree species, as they deliver multiple services like all trees, but they can also produce fruits for animals, for farmers themselves, or for diversification of the farm incomes. However, choosing fruit trees is a tougher job than choosing timber species, as fruit trees (including their rootstocks) are more sensitive to a great diversity of pests and diseases, making it much more complicated to have regular yield when an inappropriate choice is made. Additionally, few advisors are able to help the farmer willing to plant fruit trees, in regions where such crops are not well spread. For these reasons, members of the French 'RMT AgroforesterieS' (a network devoted to information-sharing about all agroforestry systems) decided to work on the development of a decision tool to help farmers in species, rootstocks and cultivars choice. The tool was developed in 2021 by a group of students and fruit experts of the RMT network. The tool has specific objectives, as it is dedicated to (i) systems where fruit is a secondary production, not the main one, and (ii) agroforestry systems where determinants will not only be technical, but also socioeconomical.

The main knowledge gathered helped to identify and order the main criteria that will have the most decisive impact on choice relevance. Soil quality, climatic conditions (drought/frost risk level) or technical complexity to grow the species have been sorted as most determining criteria. Other socio-economical traits had to be considered and added for such specific systems: labour requirement, harvest period, post-harvest constraints… as farmers have to manage other crops or livestock in priority. All criteria have been gathered and weighed under an algorithm, using R software. As a final output, experts gave a cultivar list for all 6 rosaceous species considered (apple, pear, peach, apricot, cherry, plum), according to the region where the grower is from: these cultivars are already known to be rather tolerant to pests In France, farmers are getting more and more interested in integrating trees in their agroecosystems, thanks to the French plan called « Plantons des haies » giving substantial subsidies for plantations. Breeders also consider the multiple interests of trees for broilers or even ruminants. Fruit trees are among the preferred tree species, as they deliver multiple services like all trees, but they can also produce fruits for animals, for farmers themselves, or for diversification of the farm incomes. However, choosing fruit trees is a tougher job than choosing timber species, as fruit trees (including their rootstocks) are more sensitive to a great diversity of pests and diseases, making it much more complicated to have regular yield when an inappropriate choice is made. Additionally, few advisors are able to help the farmer willing to plant fruit trees, in regions where such crops are not well spread. For these reasons, members of the French 'RMT AgroforesterieS' (a network devoted to information-sharing about all agroforestry systems) decided to work on the development of a decision tool to help farmers in species, rootstocks and cultivars choice. The tool was developed in 2021 by a group of students and fruit experts of the RMT network. The tool has specific objectives, as it is dedicated to (i) systems where fruit is a secondary production, not the main one, and (ii) agroforestry systems where determinants will not only be technical, but also socioeconomical. The main knowledge gathered helped to identify and order the main criteria that will have the most decisive impact on choice relevance. Soil quality, climatic conditions (drought/frost risk level) or technical complexity to grow the species have been sorted as most determining criteria. Other socio-economical traits had to be considered and added for such specific systems: labour requirement, harvest period, post-harvest constraints… as farmers have to manage other crops or livestock in priority. All criteria have been gathered and weighed under an algorithm, using R software. As a final output, experts gave a cultivar list for all 6 rosaceous species considered (apple, pear, peach, apricot, cherry, plum), according to the region where the grower is from: these cultivars are already known to be rather tolerant to pests


Raphaël Paut
email: raphael.paut@inrae.fr

Léa Garreau

Rodolphe Sabatier

Marc Tchamitchian

A global experimental dataset of intercropping and agroforestry studies in horticulture

Keywords: agroforestry, intercropping, horticulture, data

Intercropping fruit trees and vegetable crops in mixed systems is attracting more and more farmers in Europe (Léger et al. 2018;Lauri et al. 2019). There is also an increasing amount of research being conducted to understand the functioning of these complex systems (Do et al. 2020;Paut et al. 2021). Yet, so far, no dataset has been provided for a systematic synthesis of existing data on intercropping and agroforestry experiments in the specific field of horticulture. Therefore, the aim of the present work is to present a global experimental dataset of intercropping and agroforestry studies in the specific field of horticulture (fruit trees and vegetable crops).

The dataset includes, to date, results of field experiments published in 92 articles published between 1980 and 2017 covering 695 experiments worldwide (crop*field site* treatment combinations). The selected experiments were carried out over five continents on 11 trees and 67 vegetable crop species (Figure 1).

We plan to complete the dataset with articles published after 2017. The dataset contains measurements of variables that are considered relevant for the analysis of agroforestry and intercropping experiments: i. Basic publication information: Title, Authors, year, Country, Site (latitude, longitude).

ii. Crops and trees information: Species, Botanical family, Functional group (legume/non-legume), Plant functional type (C3/C4). CROP AND GRASSLAND PRODUCTIONS reusable and was designed to be updated. We also aim to establish a set of recommendations to design experiments on intercropping. We believe this will provide valuable data for more comprehensive analyses of intercropping and agroforestry systems in horticulture worldwide.

  iii. Parameters of the experiment: Soil texture and pH, crop management (greenhouse/field conditions), type of fertilisation (organic/inorganic), Use of pesticides, intercropping pattern (strip, row, mixed), Intercropping design (replacement, additive), Density of each species, Sowing and harvesting dates, N fertilisation, Yields of pure and mixed stands, Standard deviation of yields, Land Equivalent Ratio (LER).Data was extracted from tables, figures and text. Digitising the data from published scatter plots in the literature was performed with the WebPlotDigitizer tool (https://automeris.io/WebPlotDigitizer/). These data were assembled and harmonised according to the International System of Units. The dataset is

Figure 1 .

 1 Figure 1. Latitude and longitude coordinates of the field sites included in the database. The Köppen-Geiger climatic classification was used to link each field site to a grid size with a resolution of 0.50 degrees of latitude by 0.50 degrees of longitude. Each Köppen-Geiger climatic subzone is indicated by a colour gradient


A M Schott

M Viprey

A Termoz

J Haesebaert

M Chakir

C Della Vecchia

N Nighoghossian

S Cakmak

S Meyran

B Ducreux

C Pidoux

T Bony

M Douplat

V Potinet

A Sigal

L Derex

Evolution de la prise en charge de l'AVC à la phase aiguë en Rhône-Alpes: Etude populationnelle avant-après

Résultats:

Le nombre de patients victimes d'un IC confirmé était de 696 patients en 2007 et de 717 en 2016. Les tableaux 1 et 2 présentent les caractéristiques principales des patients et les modalités et temps de prise en charge. Après ajustement pour les facteurs de confusion (âge, sexe, régulation par le 15, période avant-après) :

-Les patients victimes d'un IC en 2007 avaient 55% moins de chances de recevoir un traitement de reperfusion que les patients victimes d'un IC en 2016 (ORaj=0,45, IC 95% 0,31-0,65).

-Le genre n'était pas significativement associé avec la probabilité d'accès à un traitement de reperfusion (ORaj=0,92, IC 95% 0.67-1.28).

-Un an d'âge supplémentaire était associé à une réduction de 1.2% de la probabilité de recevoir un traitement de reperfusion (ORaj=0,98, IC 95% 0,97-0,99) -Les patients qui n'avaient pas appelé le 15 avaient une chance d'accéder aux traitements réduite de 83% (ORaj=0,17, IC 95% 0,11-0,25),

AVC69

o

  de patients présentant une suspicion d'AVC à la phase aiguë par le SAMU 69 ou dans un des services d'urgences ou UNV du Rhône . Période avant: Cohorte AVC 69 -2007 Recrutement de 1306 patients sur une période de 7 mois (06/11/2006 et 06/06/2007) o Période après: Cohorte STROKE 69 -2016 Recrutement de 1764 patients sur une période de 7 mois (06/11/2015 au 06/06/2016) Mesures mises en oeuvre dans le Rhône: -Création d'une 2 ème UNV et augmentation des lits d'UNV de 6 à 18, -Implémentation d'un programme de formation des professionnels des SAU (étude AVCII Trial), -Campagnes locales de sensibilisation grand public (étude REACT), -Actions moins formalisées comme la formation des personnels du centre 15, des pompiers, stages immersifs en UNV de médecins urgentistes Tableau 1: Caractéristiques des patients et des modalités de prise en charge Tableau 2 : Délais de prise en charge Conclusion : Les actions entreprises dans le Rhône entre 2007 et 2016 ont amélioré l'accès des patients victimes d'un IC à un traitement de reperfusion avec une augmentation de 9% à 23%. Ceci est essentiellement du à la réduction des délais intra-hospitaliers et l'augmentation des appels au 15 mais pas à une diminution du délai entre les symptômes et l'appel au 15.

Med (Q1-Q3) min N=696 STROKE69 Med (Q1-Q3) min N=717 p Symptômes -Appel au 15

  Rationnel : En 2007 parmi tous les patients admis pour une suspicion d'AVC dans un service d'accueil des urgences (SAU) ou une unité neurovasculaire (UNV) du Rhône (cohorte AVC 69), seulement 9% des patients victimes d'un infarctus cérébral (IC) avaient eu accès à un traitement thrombolytique. En 2016, nous avons conduit une étude similaire, STROKE 69, pour évaluer l'efficacité de toutes les mesures mises en oeuvre entre 2007 et 2016 sur le taux d'accès des patients à un traitement de reperfusion (thrombolyse et/ou thrombectomie) et les délais pré et intra-hospitaliers.

		AVC 69	STROKE 69
		2007	2016	p
		N=696	N=717
	Femme (n, %)	362 (52)	360 (50)	0.50
	Age (moyenne, écart type)	75,1 (14)	74,7 (15)	0.99	Appel au 15 --Admission	30 (5-75) 67 (54-93)	41 (13-166) 64 (50-79)	0.070 0.098
	Appel au 15 (n, %)	215 (40)	488 (68)	<0.000 1	Symptômes --Thrombolyse et/ ou	196 (174-245)	155 (125-199)	<0.0001
	Admission directe en UNV Thrombolyse et/ou Thrombectomie (n, %)	57 (8) 60 (9)	214 (30) 167 (23)	<0.000 1 <0.000 1	thrombectomie Imagerie -Thrombolyse et/ ou thrombectomie Admission -Thrombolyse et/ ou thrombectomie	70 (60-98) 110 (55-166)	33 (21-48) 50 (29-78)	<0.0001 0.0002


Impact de la corrélation entre les durées et les débits des flux sur les propriétés statistiques du traffic agrégé dans les réseaux

Résumé : Depuis la découverte, en 1993, de la longue mémoire dans les traces de trafic réseaux, de nombreux modèles, basés sur les distributions à queue lourde de certaines propriétés du trafic à l'échelle flux, ont été proposés pour reproduire et interpréter cette propriété. Cependant, aucun de ces modèles ne tient compte de la corrélation existant entre les débits et les durées des flux. Dans ce travail, nous étendons les modèles existants pour inclure cette corrélation. Le modèle que nous proposons est basé sur un processus de Poisson dans le plan décrivant la structure du trafic à l'échelle flux. Il nous permet de calculer analytiquement la fonction d'auto-covariance du trafic agrégé, et de montrer que celui-ci possède la propriété de longue mémoire, avec un paramètre de Hurst différent. Dans le cas décorrélé, notre modèle est cohérent avec les modèles existants, et prédit le même paramètre de Hurst. Nous montrons également qu'une pseudo longue dépendance peut apparaître avec un autre paramètre de Hurst dans le cas de débits de flux très variables. Nos hypothèses de modélisation sont enfin favorablement confrontées à des traces de trafic web réel. Deep understanding of network-traffic properties is essential for Internet Service Providers to optimally control the traffic in order to offer users the best Quality of Service. In this context, a lot of research has been recently focusing on mathematical modeling of network traffic, especially from a statistical viewpoint. However, comprehensive modeling of all the characteristics of the traffic is a very arduous problem for it encompasses several difficulties of different natures, such as transport protocols, control mechanisms, and complexity due to the users behavior. The design of simple models, yet rich enough to account for the most important characteristics observed in the traffic, is thus a major issue for many industrial applications such as network dimensioning, optimization, control, and prediction.

Mots

A major breakthrough in this direction has been the discovery in 1993 of the self-similar nature of aggregate time series at large time scales [?, ?]. Following up, several models was proposed that posit the heavy-tailed nature of the activity periods, or ON periods, as a plausible explanation for the origin of the self-similar property; and clarify the relation between the tail exponent α ON and the Hurst parameter H:

H = 3 -α H 2 , ( 1a
)
where

α H = min(α ON , 2), or min(α ON , α OFF , 2), (1b)
depending on wether the model include some heavy-tailed OFF periods of tail index α OFF or not. Since then, the self-similarity property has been shown to have a major impact on QoS, and such models, based on a flow-level characterization of the sources, are now commonly used to generate realistic traffic traces. However, all of the models proposed up to now that lead to relation (1) entail a common simplified feature: the tail indices of the flow-size distribution α SI , and of the flow-duration distribution α ON are identical, so that both can interchangeably be used in relation (1). Even in models randomly drawing the flow rates, this equality holds due to the assumption that flow rates are independent of flow sizes (or flow durations). Such an assumption might not hold in real internet traffic, and the equality between α ON and α SI would then fall down. This have been observed since 1997 [?]. Moreover, the correlation between flow rates and durations is likely to become more prominent in the future Internet, due to the emergence of new mechanisms such as the FTTH increasing the rates spectrum, and flow-aware control procedures which might strongly correlate the achieved flow rates to flow properties such as its duration or size. In this case, relation (1) is not ensured to be reliably usable to predict the Hurst parameter of the aggregate traffic. In this context, the development of a model including the correlation between flow rates and durations is a major challenge to accurately predict the Hurst parameter of the aggregate traffic. To the best of our knowledge, no such model have been proposed yet.

In this paper, we propose an extension of existing models, which takes into account the correlation between flow rates and durations (responsible for the different tail indices α ON and α SI ). Based on a representation of the flows as a planar Poisson process, we analytically calculate the auto-covariance function of the aggregate traffic's instantaneous bandwidth and deduce the resulting Hurst parameter. We first briefly review existing models that reproduce the longrange dependence property in Section 2. Then we extend in Section 3 these models to include the correlation between flow rates and durations and show how to predict the resulting long-range dependence of the aggregate traffic. In Section 4, we experimentally demonstrate the ability of our model to accurately predict the Hurst parameter, based on a real traffic trace.

Related work

There exists a large number of models able at reproducing the long-range dependence property [?, ?]. We describe here only those which are explicitly relying on the notion of flows, particularly meaningful in network-traffic applications. Following up Mandelbrot's idea, these models rely on the introduction of a heavy-tailed distribution of infinite variance. We distinguish two categories of such models, which mainly differ in the flow arrival process: the infinite source Poisson models and the renewal models.

Infinite source Poisson models. In the simplest version of the infinite source Poisson models, flows arrive at the link as a Poisson process of rate λ, and transmit data at a fixed rate of 1 during a heavy-tailed random time of tail index α ON . It is also known as the M/G/∞ model originally considered in [?]. Then, two different limiting regimes have to be considered when studying the cumulative bandwidth at large time scales T : if λ goes to infinity first, then T , and a proper rescaling is applied, the resulting process is a fractional Brownian motion of Hurst parameter as in equation (1) (see e.g.surveys in [?, ?]). In the opposite limiting regime where the time scale T goes first to infinity, the resulting process is a Lévy motion. In the case where both λ and T goes to infinity in the same time, [?] provides conditions on the ratio between the growth rate of both quantities that ensure one, or the other limiting regime (see also [?] where an intermediate case between these two conditions is studied). There exist many variants of this simple M/G/∞ model, which all rely on the INRIA same mechanism of a Poissonian arrival of some heavy-tailed flows and mainly differ on the way data is transmitted within a flow (see also a survey of infinite Poisson models in [?], and the references in [?]). In [?], the authors consider an infinite source Poisson model with a general form of the "workload function" inside a flow and establish (for the first time) the Gaussian limit result. A similar model is considered in [?] (where the "workload function" is called "transmission schedule") and the other limiting result is shown (the Lévy motion). In [?], the Poisson shot-noise model is developed. This model is basically very similar to the two models previously mentioned. The "workload function" or "transmission schedule" is now called "shot", but "shots" still arrive according to a Poisson process. An application of this model is proposed in [?] where the shot shape is representative of the AIMD mechanism with Poisson losses. In [?], the authors consider a model where the rate within a flow is constant, but the rate of a flow is a random variable. The flow sizes are drawn at random according to a heavy-tailed distribution of tail index α SI ; and the rates are drawn at random independently of the sizes, according to any finite-mean distribution. Flow duration (flow size divided by flow rate) is then heavy-tailed with α ON = α SI and they find a long-range correlation structure of Hurst parameter H as in relation (1). A similar model where the durations and rates are drawn independently is proposed in [?]. The last model we mention in this infinite Poisson model category is the Cluster Point Processes (CPP) model proposed in [?]. In this model, the discrete approach of point processes is used. Flows are "clusters" of points and the flow arrival time is the time of the first packet of the cluster. The flow-arrival process is Poisson. The number of points in the cluster (the flow size) is heavy-tailed of tail index α SI , and points (packets) within a cluster (flows) follow a renewal process with some inter-renewal distribution determining the mean flow rate. Based on results on point processes (see [?, ?]), the authors calculate the aggregate traffic's spectrum and again find long-range dependence of Hurst parameter H as in relation (1a) where α H = α SI . As in the model of [?], the flow duration in this model is heavy-tailed of tail index α ON = α SI . As already observed in [?], this might not hold true in real Internet traffic. For example, slow-start effects might imply higher rates for larger flows. In this case, flow rate is correlated to flow size and the tail indices of the flowduration and the flow-size distributions are different. In [?], the authors suggest as future work to introduce multi-class CPP, for example with a class for small flows (mice) having small rates, and another class for large flows (elephants) having high rates. Such a work, however, has not been done and the following question, that we address in this work, remained open: What happens if flowduration and flow-size distributions have different tail indices? Which of these tail indices, if any, governs the long-range dependence of the aggregate traffic?

Renewal models. Another, closely related, class of models generating longrange dependence is the class of renewal models. They are based on the same general setting firstly introduces by Mandelbrot in 1969 [?] in an economical context. Each source, emitting only one flow at a time, is modeled as a renewal reward process, where the inter-renewal time (i.e.the interval between two consecutive flows) is heavy-tailed; and the reward is the rate of the flow. Many variants of such models have been considered, mainly differing in the distribution of the reward [?, ?, ?] (see also [?]). We present here only briefly the particular variant where the reward strictly alternate between 0 and 1: the ON/OFF model [?]. This model allows including the notion of idle time between the transmission of two flows through the OFF periods. Moreover, ONand OFF-periods distributions can have different tail indices α ON and α OFF . Then, it is shown in [?] that the same two limiting regimes as for the infinite source Poisson models yield the same limit processes. In particular, in the Gaussian limit, a fractional Brownian motion is found, whose Hurst parameter satisfies equation ( 1). The difference of this model, as compared to infinite source Poisson models is that it can account for long-range dependence induced by heavy-tailed idle times. The flow arrival process is also no longer, in general a Poisson process. However, when α ON > α OFF , which is a common case in network traffic (see e.g.[?]), the results in terms of predicted Hurst parameter are the same for both models.

Planar point process setting. In our model, we choose to use the infinite source Poisson model, and we represent the flows as a planar Poisson process. Such a setting has been introduced in the context of multifractal analysis in [?] to extend binomial multiplicative cascades. It is also used in [?] also in the context of multiplicative processes. To the best of our knowledge, this setting has been used in the context of additive processes to study long-range dependence only in the previously mentioned papers [?, ?]. ?, ?]. This means in particular that the sequence (D i ) i≥0 is independent. The measure Λ determines the mean number of points in a part of the plane. For example, Λ(C(t 1 )) is the mean number of points in the cone C(t 1 ) (see Figure 1). In the case where the measure has a constant density (or intensity), the point process is called homogeneous. Here, the measure contains the information of the flow duration distribution. To account for a heavy-tailed flow-size distribution, we use the following form for the measure Λ of an elementary square of size (dt, dd) centered on (t, d):

A model accounting for the correlation between flow rates and durations

Definitions and notations

Λ(dt, dd) = Cdtdd d αON+1
, where C > 0, and α ON > 1.

(2)

Due to the dependence of Λ in the variable d, the point process is non-homogeous (the density of points is higher for smaller values of d). However, due to the independence of Λ in the variable t, the arrival time of a flow T i is independent of its duration D i . In our case where the x-axis represents the time, this specific form of the measure Λ ensures the stationarity of the resulting traffic. To avoid integrability problems for d around zero, we also set a minimal flow duration ε, under which the distribution is 0. This threshold will not play any role in our study because we concentrate on long-range dependence. A similar setting is used for example in [?], where the behavior of the measure around d = 0 has a great role because the authors focus on small-scale multifractality (note however, that in this paper the authors study a multiplicative process which differs from our additive process considered here). Each flow emits data at a constant rate R i during time D i . The rates are random variables such that sequence (R i ) i≥0 is independent. However, contrarily to previously proposed models (see in particular [?, ?]), we do not assume that R i is independent from D i . Instead, we permit some correlation between the flow rates and durations that we will precise later in this section, when it becomes necessary to pursue the calculations. The flow size, that we denote S i is then determined by

INRIA λ W = E{W t 1 ,t 2 } D i T i t 1 t 2 R i time t duration d W t 1 ,t 1 µ W = E{W t 1 ,t 1 } W t 2 ,t 2 µ W = E{W t 2 ,t 2 } C(t 1 ) C(t 2 ) W t 1 ,t 2
S i = D i R i .
We consider two time instants t 1 < t 2 and introduce the following notations (see Figure 1):

W t1,t1 = (Ti,Di)∈C(t1)\C(t2) R i , (3)
W t1,t2 = (Ti,Di)∈C(t1)∩C(t2) R i , (4)
W t1 = W t1,t1 + W t1,t2 = (Ti,Di)∈C(t1) R i . (5)
The variables W t2,t2 and W t2 are defined similarly. The random variable W t1 is directly the instantaneous throughput at time t 1 (the sum of the rates of the flows active at time t 1 ). The random variables W t1,t1 and W t1,t2 are intermediate random variables, useful for the calculations developed in the following. We can interpret these random variables as follows: W t1,t1 corresponds to the traffic of the flows that are active at time t 1 but not any more at time t 2 , whereas W t1,t2 corresponds to the traffic of the flows active at time t 1 and still active at time t 2 . For clarity purposes, we introduce the following notations for the mean:

λ W = E{W t1,t2 } = Λ(C(t 1 ) ∩ C(t 2 )), (6)
µ W = E{W t1,t1 } = E{W t2,t2 } = Λ(C(t 1 )\C(t 2 )) = Λ(C(t 2 )\C(t 1 )), (7)
where simple calculations show that the last equality Λ(C(t 1 )\C(t 2 )) = Λ(C(t 2 )\C(t 1 )) holds due to the uniform distribution of the measure Λ in time. Finally, note that with these notations, we have:

E{W t1 } = E{W t2 } = λ + µ. ( 8
)
The setting described here allows us to compute the autocovariance of the instantaneous bandwidth W t , in order to evaluate the Hurst parameter, as we shall see now (recall that a stationary process is called long-range dependent of Hurst parameter H if its auto-covariance function decrease like τ 2-2H when τ goes to infinity, where 1/2 < H < 1).

Computation of the instantaneous bandwidth correlation

We are interested in the computation of the autocovariance function of the instantaneous bandwidth, that is

E{W t1 W t2 } -E{W t1 }E{W t2 }. Proposition 3.1. If the planar point process {(T i , D i )} is Poisson, then E{W t1 W t2 } -E{W t1 }E{W t2 } = E{W 2 t1,t2 } -E{W t1,t2 } 2 = Var{W t1,t2 } (9) Proof.
Using the notations introduced above, we have:

E{W t1 W t2 } = E{W t1,t1 W t2,t2 }+E{W t1,t1 W t1,t2 }+E{W t1,t2 W t2,t2 }+E{W 2
t1,t2 }. Due to the Poisson assumption, the random variables W t1,t1 , W t2,t2 , W t1,t2 are mutually independent, so that

E{W t1 W t2 } = E{W t1,t1 }E{W t2,t2 } + E{W t1,t1 }E{W t1,t2 } + E{W t1,t2 }E{W t2,t2 } + E{W 2 t1,t2 } = µ 2 W + 2λ W µ W + E{W 2 t1,t2 } = (λ W + µ W ) 2 -λ 2 W + E{W 2 t1,t2 } INRIA
directly giving equation ( 9) in view of the definition of λ W (equation ( 6)) and of equation ( 8).

Proposition 3.1 shows that the autocovariance depends only on the variance of the traffic due to the flows in the intersection of the cones C(t 1 ) and C(t 2 ). This had been noticed already in [?, ?], though the authors were focusing on small-scale properties. To complete the computation of the autocovariance function, we then need to compute Var{W t1,t2 }.

If the flow rates were constant and equal to 1, then W t1,t2 would simply be the number of points in C(t 1 ) ∩ C(t 2 ). Since the point process is Poisson, the variance Var{W t1,t2 } would then simply be λ W and the autocovariance would be determined by the value of λ W (recall that the variance of the count process associated with a Poisson process is equal to its mean). Before proceeding with the calculations in a more general case, we introduce additional useful notations. We denote by N the random variable corresponding to the number of points in C(t 1 ) ∩ C(t 2 ). We denote its mean by

λ N = E{N }. ( 10
)
Note that in the case where the rate is always 1, we simply have λ W = λ N . The next proposition gives the general form of the autocovariance function, without specifying the correlation between R i and D i yet.

Proposition 3.2.

E{W t1 W t2 } -E{W t1 }E{W t2 } = λ N E{R 2 i }, ( 11
)
where it is implicitly understood that the expectation is computed in C(t 1 )∩C(t 2 ).

Proof. To evaluate the value of Var{W t1,t2 }, we successively compute the values of E{W t1,t2 } and E{W 2 t1,t2 }.

E{W t1,t2 } = E{E{W t1,t2 |N }} = ∞ k=1 E{ k i=1 R i |N = k}P(N = k) = ∞ k=1 k i=1 E{R i |N = k}P(N = k) = ∞ k=1 kE{R i }P(N = k) = λ N E{R i }. E{W 2 t1,t2 } = E{E{W 2 t1,t2 |N }} = ∞ k=1 E{( k i=1 R i ) 2 |N = k}P(N = k) By independence of the sequence (R i ) i E{( k i=1 R i ) 2 |N = k} = E{( k i=1 R i ) 2 } = kE{R 2 i } + (k 2 -k)E{R i } 2 , so that E{W 2 t1,t2 } = ∞ k=1 (kE{R 2 i } + (k 2 -k)E{R i } 2 )P(N = k) = ∞ k=1 kE{R 2 i }P(N = k) + ∞ k=1 k 2 E{R i } 2 P(N = k) - ∞ k=1 kE{R i } 2 P(N = k) = λ N E{R 2 i } + (λ N + λ 2 N )E{R i } 2 -λ N E{R i } 2 = λ N E{R 2 i } + λ 2 N E{R i } 2
Until that point, we still have not used the precise form of the measure Λ (equation ( 2)), and the precise form of the correlation between R i and D i . The result of Proposition 3.2 depends only on the Poisson and independence of the sequence (R i ) i assumptions. It shows that the autocovariance function is the product of two terms: λ N , the mean number of points in C(t 1 ) ∩ C(t 2 ), and E{R 2 i }. The first term (λ N ) depends only on the measure Λ and is easily obtained via a simple integration (Proposition 3.3). To compute the second term (E{R 2 i }), we need in addition to precise the correlation between R i and D i (Proposition 3.4).

Proposition 3.3. If the measure Λ has the form of equation ( 2), then

λ N = C 1 α ON (α ON -1) (t 2 -t 1 ) -αON+1 (12)
Proof.

λ N = ∞ d=t2-t1 t1 t=t1-(d-(t2-t1)) Λ(du, dv) = C ∞ d=t2-t1 (d -(t 2 -t 1 )) 1 d αON+1 dd = C ∞ d=t2-t1 1 d αON dd -C(t 2 -t 1 ) ∞ d=t2-t1 1 d αON+1 dd = C 1 α ON -1 (t 2 -t 1 ) -αON+1 -C 1 α ON (t 2 -t 1 ) -αON+1 = C 1 α ON (α ON -1) (t 2 -t 1 ) -αON+1
Finally, we now specify the correlation between R i and D i to calculate the term E{R 2 i }. Our goal here is to specify a correlation which indeed leads to the different tail indices of the flow-size and flow-duration distributions. We have already mentioned that taking R i as a random variable independent of D i would lead to the same tail indices, independently of the distribution of R i , provided that its mean is finite (see for instance [?, ?]). The different tail indices α SI and α ON can then only come from correlation between R i and D i . The simplest choice would be to deterministically take for each flow: R i = D β-1 i ,

INRIA

where β = αON αSI . In this case, we would have S i = D β i or equivalently

D i = S 1/β i
. This would effectively lead to the different tail indices α ON and α SI for the flow duration and size distributions. However, this assumption of a deterministic rate for a flow of a given duration is not realistic. Instead we choose a model where the conditional expectation and variance of the rate given the duration are specified. This model is given in the next proposition.

Proposition 3.4. Assume that E{R i |D i } = KD β-1 i
, where K is a constant and Var{R i |D i } = V , where V is a constant. We denote

α ′ = α ON -2(β -1). ( 13
) If α ′ > 1, then E{R 2 i } = 1 λ N CK 2 1 α ′ (α ′ -1) (t 2 -t 1 ) -α ′ +1 + V. ( 14
)
Proof.

E{R 2 i } = E{E{R 2 i |D i }} = K 2 E{D 2(β-1) i } + E{Var{R i |D i }}
The same kind of integration as for the proof of Proposition 12 give

E{D 2(β-1) i } = 1 λ N C 1 α ′ (α ′ -1) (t 2 -t 1 ) -α ′ +1
, while we clearly have E{Var{R i |D i }} = V , completing the proof.

We are now able to state the final proposition establishing the decrease of the autocovariance function with t 2 -t 1 , and then the long-range dependence of the process (W t ) t .

Proposition 3.5 (Autocorrelation function and long-range dependence of the process

(W t ) t ). If E{R i |D i } = KD β-1 i and Var{R i |D i } = V , where K, V are constants and α ′ = α ON -2(β -1) > 1, then E{W t1 W t2 } -E{W t1 }E{W t2 } = CK 2 1 α ′ (α ′ -1) (t 2 -t 1 ) -α ′ +1 +CV 1 αON(αON-1) (t 2 -t 1 ) -αON+1 . (15)
The process (W t ) t is then (asymptotically) long-range dependent with Hurst parameter H following equation ( 1a) (H = 3-αH 2 ), where α H = min(α ′ , α ON ).

(16)

Proof. The first part is immediate from Propositions 3.2, 3.3 and 3.4. For the Hurst parameter, recall that a process is long-range dependent of Hurst parameter H if its autocovariance function algebraically decreases like (t 2t 1 ) 2H-2 .

Proposition 3.5 is our main result for this section. It establishes the longrange dependence of the instantaneous bandwidth and gives the Hurst parameter. Before going back to the in2p3 trace and verifying the consistence of our model choice with the data, we make a few general remarks on the result of Proposition 3.5, and the choice of the model.

Let us first comment on specific values of β = αON αSI .

If β = 1 (α SI = α ON ): This case corresponds to the classical case of Taqqu's model and other infinite source Poisson models, where the rate is not correlated to the duration and we have the same tail indices for the flow-size and flow-duration distributions. In this case, the result of Proposition 3.5 reduces to Taqqu's relation (1a), where the tail index controlling the long-range dependence α H is the unique tail index α SI = α ON (bear in mind that there are no heavy-tailed OFF times here because we used a Poisson flow arrival).

If β > 1 (α SI < α ON ): This is the case of the in2p3 trace, where the rate increases in average with the duration of the flows. In this case the tail index controlling the long-range dependence is α H = α ′ < α ON . The long-term correlations are then stronger than in the case of constant rates. Note that, depending on the value of β, we can have either α H < α SI , α H > α SI or α H = α SI . It means that the tail index controlling the long-range dependence does not necessarily lie between α SI and α ON , but can be smaller than both these tail indices.

If β < 1 (α ON < α SI ): This case corresponds to a situation where the rate decreases in average with the duration of the flows. This could happen for instance with some scheduling policy giving priority to small flows. In this case the tail index controlling the long-range dependence is α H = α ON < α ′ .

In the case where V = 0 (this corresponds to the deterministic case

R i = D β-1 i
), the tail index controlling the long-range dependence is always α H = α ′ , even if β < 1. Note that depending on the values of β, we can have in this case surprising situations similar to the situations discussed for β > 1 where the tail index controlling the long-range dependence α H is greater than the two tail indices α SI and α ON . Note that the calculation of the autocovariance function of equation ( 15) could easily be performed with some other (non-constant) forms of the conditional variance Var{R i |D i }. For example, an algebraically increasing variance (Var{R i |D i } = D γ i ) would lead to another index α ′′ = α ON -γ, possibly controlling the long-range dependence for some set of parameters.

The long-range dependence that we have discussed until now is asymptotic long-range dependence. However, depending on the values of the constants K and V , we can observe situations where within some intermediate range of scale, the term of the autocovariance function (15) with the larger exponent dominates. For example, in the case β > 1, the first term of the autocovariance asymptotically dominates, but if the value of V is very large and the value of K is not too large, the second term can dominate within some scale range, thus leading to pseudo long-range dependence controlled by the tail index α ON . Note that the two constants K and V do not have the same units and cannot be compared directly to one another. Instead, the whole terms of equation ( 15) have to be compared.

Let us finally mention that to "specify enough" the correlation between the flow rates and durations, in order to compute the autocovariance function of equation ( 15), we had to specify only the first two conditional moments E{R i |D i } and Var{R i |D i }. This is not surprising since the autocovariance function is a second-order moment quantity, and we found only second-order self-similarity (i.e., long-range dependence). To investigate finer statistical properties of the instantaneous bandwidth, we might need to specify the entire conditional probability P(R i = r|D i = d).

INRIA 4 Results and discussion

Trace description

In this section, we use a real traffic trace to validate the ability of our model to predict the aggregate-traffic Hurst parameter. The trace is acquired at the output link of the in2p3 research center (Lyon, France), with the capture tool MetroFlux [?]. The traffic is captured from the VLAN corresponding to RE-NATER 1 web traffic. This VLAN is encapsulated in the 10 Gbps output link of in2p3. Although we captured more than one day of traffic, we restrict our trace to a 30 minutes stationary trace, corresponding to the incoming traffic between 3pm and 3:30pm on January 18, 2009. The mean throughput in this period is 127.3 Mbps. This traffic mainly corresponds to web traffic, but has the particularity of including a larger number of elephants than usual traffic. Since the in2p3 is a nuclear-physics research center, these large transfers are likely to correspond to the transfer of experiment results from experiment centers like CERN.

Figure 2 shows the flow-size and flow-duration distributions. While both distributions clearly appear to be heavy-tailed, they exhibit sharply different tail indices: α SI = 0.8544 and α ON = 1.1994 respectively (the tail indices are estimated with the recent wavelet method proposed in [?]). Several explanations could be posited to interpret why longer flows achieve higher rates, thus explaining these different tail indices. A plausible explanation might lie in transient effects at the beginning of each flow. Another explanation could lie in the variable locations of the downloaded files, if for example users tend to download large flows from closer locations (thus achieving higher rates more quickly because of smaller RTTs). The largest files might also be intentionally stored on servers with the highest capacities. These are only hypotheses and we do not elaborate further here on the origin of the different tail indices for the flow-size and flow-duration distributions, which is out of our scope.

To investigate the long-range dependence of the aggregate traffic, we use the wavelet method described in [?]. This method provides a robust estimation of the Hurst parameter, with an estimation of the confidence interval in the Gaussian case. Here, we observed a fairly Gaussian traffic, with a kurtosis of 3.15 at scale 10 ms. Figure 3 displays the log diagram of the aggregate traffic bandwidth in time windows of size ∆ = 10 ms. We clearly observe a linear behavior in the coarse scales, characteristic of a scaling behavior, with an estimated Hurst parameter of H = 0.901 ± 0.045. The estimation is performed in the scale range [0.64 -40.96] s (lying far beyond the mean ON time of 0.12 s) but the same scaling behavior appears to extend to the scale 100 s. This Hurst parameter shows a good agreement with relation (1) of Taqqu's ON/OFF model, where α H = α ON would be the tail index governing the long-range dependence. However, as already mentioned, the trace studied here is not well modeled by the ON/OFF model with constant rates, and we now show that the model developed in previous section match well the real trace and allows us to correctly predict this Hurst parameter. Note finally that we observed on the trace a coarsescale parameter of 0.65 for the flow arrival process. Such a value is unlikely to explain the long-range dependence of the aggregate traffic. Although not strictly satisfied, the Poisson flow arrival assumption made in our model should

Confrontation of the model with the trace

Our goal here is to confront the model we proposed in previous section to the real data, especially regarding the assumptions we made on the correlation between flow rates and durations in Proposition 3.5. From the in2p3 trace, we have Figure 4 shows the conditional moments E{R i |D i } and Var{R i |D i }, as functions of D i , estimated from the in2p3 trace. We plotted lines corresponding to the model of Proposition 3.5: E{R i |D i } = KD β-1 i and Var{R i |D i } = V , where K, V are constants, with the estimated value of β. These lines are vertically adjusted and show a good match in the range [1-100] s, with the values K = 10 5.2 , V = 10 12.4 . This confirms the relevance of our model in this scale range, which roughly corresponds to the scale range in which we estimated the Hurst parameter in Section 4.1. With these values of K and V , we can compare the two terms of the autocovariance function of equation ( 15). We find that in the scale range of estimation of H ([0.64 -40.96] s), the first term of equation ( 15) is always five times smaller than the second term. We conclude that the long-range dependence observed on Figure 3 is pseudo long-range dependence controlled by the second term of equation ( 15), i.e., by the tail index of the flow-duration distribution α ON . The estimated value of the Hurst parameter shows a very good agreement with relation (1a) with α H = α ON , confirming the accuracy of the model in Proposition 3.5 in predicting the Hurst parameter's value provided that the terms of equation ( 15) are appropriately compared in the scale range of interest.

INRIA E{R i |D i } 0.1s 1s

Conclusion

In this work, we developed a model taking into account the correlation existing between flow rates and durations responsible for the different tail indices α SI and α ON . Based only on specifications of the first two order moments of the conditional probability of the rates given the duration, we showed that the instantaneous bandwidth exhibits long-range dependence with a Hurst parameter as in relation (1a), where the controlling tail index α H can be either α ON or a combination of α ON and α SI , whichever is the smaller. We also showed that, in certain circumstances, pseudo long-range dependence with a different Hurst parameter can be observed in a finite scale range. In the case where the correlation between flow rates and durations vanishes, our results coincide with the usual results of classical models.

Finally, we validated the ability of our model to predict the correct Hurst parameter on a real web traffic trace. However, the non-stationarity inherent to real Internet data at large time scales (around one hour) constrained us to restrain our study to 30 minutes of traces. Consequently, we observed pseudo long-range dependence, and longer stationary traces would be required to illustrate the ability of our model to predict asymptotic long-range dependence, which, as we saw, can exhibit a different Hurst parameter.

The model proposed in this work is the first model, to the best of our knowledge, that includes the correlation between flow rates and durations. This correlation is very important. It has been observed for over a decade [?], and is likely to become even more important with the appearance of FTTH and flowaware approaches. Further theoretical developments would be needed to give a more rigorous mathematical statement of our results, in particular concerning the Gaussian limit, and a potential other limit that could arise similarly to classical models. Also, we assumed a Poisson flow arrival, which is a useful simplification to conduct the calculations, but might not always hold. Considering less restrictive flow-arrival processes would also be an interesting improvement of our results, be it only to clarify in which situations a correlated flow-arrival process can affect the aggregate traffic self-similarity.

4

  -clés : Trafic réseau, Dépendence à longue portée, Distributions à queue lourdeand notations . . . . . . . . . . . . . . . . . . . . . . 6 3.2 Computation of the instantaneous bandwidth correlation . . . . . 8 Results and discussion 13 4.1 Trace description . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 4.2 Confrontation of the model with the trace . . . . . . . . . . . .

Figure 1

 1 Figure1depicts our setting. The point process (T i , D i ) (representing the arrival time and the duration of the flows) is assumed to be a planar Poisson process of measure Λ [?, ?]. This means in particular that the sequence (D i ) i≥0 is independent. The measure Λ determines the mean number of points in a part of the plane. For example, Λ(C(t 1 )) is the mean number of points in the cone C(t 1 ) (see Figure1). In the case where the measure has a constant density (or intensity), the point process is called homogeneous. Here, the measure contains the information of the flow duration distribution. To account for a heavy-tailed flow-size distribution, we use the following form for the measure Λ of an elementary square of size (dt, dd) centered on (t, d):

Figure 1 :

 1 Figure 1: Setting of the model. Each point represent a flow. The x-coordinate represents the start time T i and the y-coordinate the duration D i . R i is the rate of the flow. At time t 1 , the active flows are those lying in cone C(t 1 ) (the left border of the cones have a slope of -1).

1Figure 3 :

 3 Figure 2: (a) Flow-size distribution -(b) Flow-duration distribution. The tail indices estimated with the wavelet method are respectively: α SI = 0.8544 and α ON = 1.1994. Straight lines materialize these slopes and have been vertically adjusted to the data. The mean ON time is 0.12 s.

Figure 4 :

 4 Figure 4: First two moments of the conditional probability of the flow rates given the flow durations. The red dots are empirical histograms estimated on the in2p3 trace. The black lines are not estimated from the empirical histograms. They have slope (left) β -1 = 0.4037 where β is estimated from the tail indices α SI and α ON -(right) slope 0. They have been vertically adjusted to fit the data in the range [1 -100] s

RR n°7100


Carlos F Crispim-Junior
email: carlos.crispim-junior@liris.cnrs.fr

Romain Guesdon

Christophe Jallais

Georges Beurier

Xuguang Wang

Florent Laroche

Stéphanie Souche-Lecorvec

Laure Tougne Rodet

A glimpse at the first results of the AutoBehave project: a multidisciplinary approach to evaluate the usage of our travel time in self-driving cars

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

A glimpse at the first results of the AutoBehave project: a multidisciplinary approach to evaluate the usage of our travel time in self-driving cars Carlos F. Crispim-Junior 1* , Romain Guesdon

Introduction:

The advances in self-driving technologies will certainly revolutionize our lives. The first robot taxi services are already present in certain cities of the world (e.g., WAYMO in San Francisco). Without surprise, most work related to self-driving cars (SDC) focuses on improving the way an SDC navigates and interacts with other road actors (e.g., pedestrians, cyclists, other vehicles). However, what happens in the cockpit of SDCs is little explored. Existing work focuses on situations where the driver is asked to take back control of the vehicle. We currently lack tools, and datasets to analyze the behaviors of the occupants of SDC in outdoor conditions, and beyond driving-related actions [START_REF] Ortega | Dmd: A large-scale multi-modal driver monitoring dataset for attention and alertness analysis[END_REF]. For instance, most computer vision datasets devoted to monitoring the cockpit of SDCs of SAE1 level 4 (high automation) consist of synthetic data [START_REF] Da Cruz | SVIRO: Synthetic Vehicle Interior Rear Seat Occupancy Dataset and Benchmark[END_REF][START_REF] Katrolia | TICaM: A Time-of-flight In-car Cabin Monitoring Dataset[END_REF][START_REF] Guesdon | Synthetic Driver Image Generation for Human Pose-Related Tasks[END_REF] or are based on laboratory recordings [START_REF] Martin | Drive&Act: A Multi-Modal Dataset for Fine-Grained Driver Behavior Recognition in Autonomous Vehicles[END_REF]. Besides, studies in human and social sciences, such as cognitive sciences and economy of transport are usually focused on declarative surveys and/or evaluated on dedicated infrastructures, such as cockpit simulators [START_REF] Malokin | How do activities conducted while commuting influence mode choice? Using revealed preference models to inform public transportation advantage and autonomous vehicle scenarios[END_REF][START_REF] Harb | What do we (Not) know about our future with automated vehicles?[END_REF].

Contribution: We propose a multidisciplinary, multi-sensor framework named Auto-Exp to evaluate how the occupants of SDCs use their travel time in real-world conditions, in particular concerning non-driving actions.

Materials and methods: The proposed framework is composed of an experimentation scenario (Fig. 1C), and an acquisition module (Fig. 1A). The experimentation scenario is composed of two parts: the free and the guided scenarios. The free scenario seeks to capture the natural (non-guided) actions of the occupants of a vehicle during their usage (or discovery in our case) of the SDC. The guided scenario takes place after the free scenario, and it consists of the realization of a series of actions using the objects put at the disposal of the participants. Its main goal is to collect a video dataset for the development of video analysis tools. The acquisition module is composed of two RGB-D sensors (Intel RealSense D435) and a GoPro camera. RGB-D sensors are attached to the top corners of the windshield and capture a lateral view of the body of people in the front seats. The GoPro camera records a close view of the upper body of the person in the driver's seat. We also ask participants to fill out declarative surveys about their internal states (IT) and their acceptance of this new means of transportation (ASDC). We applied the surveys before (ASDC+IT), during (IT only), and after the use of the vehicle (ASDC+IT).

Experiments: To evaluate the proposed framework, we carried out a four-day long experiment in July 2021 using a Renault Zoe car (electric supermini urban model) on the campus of the Ecole Centrale de Nantes (ECN) in France. The vehicle was robotized by the LS2N-ARMEN laboratory to behave as an SDC of SAE level 4 [START_REF] Dominguez | Longitudinal Dynamics Model Identification of an Electric Car Based on Real Response Approximation[END_REF]. Participants were recruited based on their declared interest to test the SDC technology during their participation in a prior, online survey on people's attitude towards SDCs [START_REF] Corvec | Autonomous vehicle: attitude, activities, and feelings[END_REF], or in response to recruitment advertisements posted to the employees and the students at the university. Participants occupied the driver's seat of the vehicle. A safety driver was present in the passenger front seat using a dedicated interface to take back control of the car in case of need. An experimenter was also present in the back seat to guide the experiment. The vehicle traveled over a predefined path over the streets and the parking of the ECN campus (Fig. 1B) and automatically adapted its trajectory in response to the road conditions. Each participant's experiment consisted of five laps in the predefined path (8 km in total). The free scenario took place on the first three laps (4.8 km on average), and the guided scenario took place during the last two laps (3.2 km on average). During the free scenario, participants were instructed to behave as they feel like. They could perform actions with objects they brought with them or with objects put at their disposal, or simply watch the experiment unfold idly. During the guided scenario, the experimenter asked participants to realize a set of actions. The instructions were brief, so the participants would have the liberty to decide how to perform the actions. We created the list of target actions based on previous studies about autonomous vehicles. [START_REF] Martin | Drive&Act: A Multi-Modal Dataset for Fine-Grained Driver Behavior Recognition in Autonomous Vehicles[END_REF][START_REF] Corvec | Autonomous vehicle: attitude, activities, and feelings[END_REF]. Examples of action instructions are "use the phone", "read a newspaper", or "listen to music". The experimenter chose randomly the order of guided actions, as well as their duration.

Results: We acquired a multidisciplinary, multi-sensor dataset composed of 29 partic-ipants (18 male/11 female) conducting non-driving actions related to a common commuting routine on an SDC of SAE level 4. Video recordings are 1 hour long and are associated with surveys about participants' internal states and attitudes toward SDC. The acquired dataset successfully depicts a large variety of sources of noise, both in terms of sensor data and action realization. The observed types of noise are common in real-world scenes but are currently absent in previous datasets. Recruited participants presented a variety of body sizes, ages, and education levels. The acquired dataset, named AutoBehave, will be released progressively on the project website2 under the form of computer vision tasks (e.g., pose estimation, action recognition) and multidisciplinary studies (e.g., the analysis of the internal states of SDC's users to estimate their comfort, the study of the value of the user's time in this new type of vehicle).

Conclusion:

We proposed a novel framework to study how occupants of SDCs use their travel time by analyzing their actions, internal states, and attitudes toward this kind of vehicle. We demonstrated its usage by collecting a multidisciplinary dataset depicting the non-driving actions of people during the test of an SDC traveling in an outdoor environment. Future work will focus on developing methods to analyze the acquired data automatically and on the study of the acquired dataset in a multidisciplinary fashion.

Figure 1 :

 1 Figure 1: Proposed multidisciplinary, multi-sensor framework. A) Acquisition module, B) Experimental path on the campus of Ecole Central de Nantes, C) Timeline and scenarios of the experimental framework.

  1 , Christophe Jallais 3 , Georges Beurier 2 , Xuguang Wang 2 , Florent Laroche 4 , Stephanie Souche-Le Corvec 4 ,Laure Tougne Rodet 1 1 Univ Lyon, Lyon 2, LIRIS UMR 5205, F-69676, Lyon, France

2 Univ Eiffel, Univ Lyon 1, LBMC UMR T 9406, F-69622 Lyon, France 3 Univ Gustave Eiffel, Univ Lyon, TS2-LESCOT, F-69675 Lyon, France 4 Univ Lyon, Université Lyon 2, LAET, F-69007, Lyon, France September 15, 2023

Society of Automotive Engineers

AURA AutoBehave: https://perso.liris.cnrs.fr/carlos.crispim-junior/aura autobehave2019.html

Acknowledgment

We developed this work in the context of the AURA AutoBehave project 2019, a French project devoted to developing novel methods to evaluate the acceptability of autonomous vehicles using a multidisciplinary approach. We thank the Region Auvergne-Rhône-Alpes for funding this research.


DEVELOPING A CONTEXTUAL DATABASE USING CENSUS AND OTHER ADMINISTRATIVE DATA The "Trajectoires et Origines 2" Contextual Database

Keywords:

BACKGROUND

• Demographic behavior is influenced by the surrounding context

• Information about the context in which people live is increasingly important in understanding demographic phenomena (Caporali et al., 2016) • National -sometimes regional-macro indicators are commonly used in demography

• Neighborhood and municipality level indicators are difficult to obtain and to access • There is a need for contextual data encompassing different geographic levels

•

  Family, health and education policies are often decided at the regional level • Policies regarding childcare facilities are often left to municipalities OBJECTIVES • A contextual database for the "Trajectories and Origines 2" (TeO2) survey, INED & INSEE (French Statistical Office) 2019-2020 • Build a set of indicators characterizing the place of residence of respondents at the time of the survey • Develop a method to build a database of geocoded contextual indicators derived from the French population census and other French administrative sources VISUALIZATION • We are developing a Shiny App for dynamic data visualization (with Arno Muller, INED)

DATA AVAILABILITY

• Make it available for record linkage to other geocoded surveys

RETROSPECTIVE & PROSPECTIVE CDB

• The current version of the CDB describes the French context in 2018

• It characterizes the survey respondents' current place of residence

• "Which context is it that matter? At which point at individual life stage?" (Goldscheider, 2000)

• Rosa Weber

Giulia Ferrari

Institut National d'Études Démographiques (INED) giulia.ferrari@ined.fr


M Loiseau

S Chataigner

R Créac'hcadec

E Lepretre

J P Court

M O Quéméré

Modelling of the creep behaviour of adhesive joints using Arcan tests

Introduction

Structural adhesive bonding is being increasingly used in civil engineering applications. This joining technique has in fact several advantages: weight reduction, facility to joint different materials, uniform stress distribution, no heat production. This last feature is particularly interesting for offshore activities, where there may be some need to connect new elements to existing structures without constraining their service properties [START_REF] Adams | Adhesive Bonding: Science, Technology and Applications[END_REF]- [START_REF] Da Silva | Modeling of Adhesively Bonded Joints[END_REF]. However, there are still some issues preventing a wider development of this technique such as the durability of structural adhesive bonding. The firm Cold Pad developed adhesively bonded fastener solutions (FIGURE 1 -(a)) whose aim is to replace welding by structural bonding, for repair or reinforcement of steel structures in offshore area. A specific installation tool (FIGURE 1 -(b)) was also developed to ensure precise control of the bonding operation between the two substrates [START_REF] Court | Bonded assembly and bonding method[END_REF]. [START_REF] Court | Bonded assembly and bonding method[END_REF] Regarding the durability issue raised by such a solution, as specific protections against humidity diffusion have been provided, creep appears to be one of the main phenomena that must be studied for the considered case. To do so, experimental investigations were carried out at real scale on the assembly in creep for two geometries of a bonded fastener, at different stress levels and under different load situations (tension and shear) using a methacrylate adhesive [START_REF] Djeumen | Creep investigations on adhesively bonded fasteners developed for offshore steel structures[END_REF]. However, these investigations remain specific to the studied fastener geometry preventing their use for creep prediction of alternative geometries. Most of existing literature focused on the creep characterisation on bulk samples of the adhesive. Yet, as interphase may influence this behaviour, studies at the scale of the adhesive joint are required.

FIGURE 1 -(a) Bonded fastener solution (b) Installation tool developed by Cold Pad

Arcan test (FIGURE 2 -(a)) presents several advantages compared to other tests at the scale of the adhesive joint: beaks closed to the bonded surfaces to reduce the edge effects (FIGURE 2 -(b)); different type of loads can be applied to the adhesive layer. It was therefore chosen for this study. The selected experimental set-up is the Arcan test TCS developed by Créac'hcadec et al. [START_REF] Créac'hcadec | A stress concentration-free bonded arcan tensile compression shear test specimen for the evaluation of adhesive mechanical response[END_REF]. The present study will focus on the analysis of the creep behaviour of bonded assemblies using the Adekit A310-1 methacrylate adhesive.

Materials and methods

The thickness of the adhesive joints is defined by the geometry of the samples. In this works, the tested adhesive thickness is 200 µm. The bonded surface is 9.5425.4 mm².To investigate the creep behaviour using the Arcan test, an adapted creep testing machine was developed (FIGURE 3). The load applied to the specimen is created thanks to gravity. It is amplified with two level arms. To control the loading speed, the upper arm is supported by a threaded rod. The rod is removed progressively to apply the load on the sample.

FIGURE 3 -Creep mechanical frame

A Digital Image Correlation (DIC) technique is often used to extract displacements data with Arcan tests [START_REF] Valès | An experimental method dedicated to the dynamic characterization of structural adhesives under drop weight conditions[END_REF]. However, for this study, Linear Variable Differential Transformer (LVDT) sensors were chosen to measure the displacements of the bonded assembly. The relevant displacements were in the normal and orthogonal directions of the bonded surfaces. A specific system was designed to fix the sensors on the beaks of the sample (FIGURE 4Erreur ! Source du renvoi introuvable.). For tensile tests, the relevant displacement is the one that is normal to the bonded surface. To validate the displacement measurements with the LVDT sensors, some quasi-static investigations up to failure were done with DIC to compare the results of displacement data with this method to the measure with LVDT sensors. These tests enabled to validate the measurement method with the LVDT sensors.

FIGURE 4 -Support for the LVDT sensors adapted to the sample

Creep test investigations

Preliminary quasi-static investigations up to failure were carried out to determine the monotonic ultimate capacities of the adhesive joints. The test was repeated at least five times. The temperature of the test was 20°C. The average ultimate capacity in tension, noted T ult , is 17.9 MPa with a standard deviation of 10%.

The creep tests in tension were performed between 50% and 70% of T ult at 20°C. FIGURE 5 and FIGURE 6 detail the results of these tests. FIGURE 5 highlights that the logarithm of the time to failure (t f ) is linearly linked to the normal stress applied to the adhesive joint. This observation was already made in other studies [START_REF] Sokairge | Long-term creep behavior of basalt fiber reinforced polymer bars[END_REF]. FIGURE 6 shows the displacements normal to the bonded surface during the creep tests as function of the normalised time. This time is defined as the ratio of the time and t f . It can be noticed that all the tests present the three creep stages. This definition of normalised time has already been used for previous creep studies [START_REF] Loiseau | Inverse identification method of adhesive creep properties from real scale investigations on bonded fastener[END_REF], [START_REF] Loiseau | Determination of Creep Behaviour of Adhesively Bonded Assembly -Application to Adhesively Bonded Steel Fasteners[END_REF].

FIGURE 5 -Time to failure as function of the normal stress for creep tests in tension FIGURE 6 -Normal displacement depending on time during tensile creep tests at 50%, 60% and 70% of T ult

To investigate the influence of temperature on the creep behaviour of the adhesive joint, one of the frames was placed in a thermal chamber. The creep tests at high temperature were performed at 35°C and 50°C, at 50% of T ult for the tensile tests. As expected, the time to failure declines sharply when the temperature increases (FIGURE 7). The logarithm of the time to failure seems linearly linked to the temperature. FIGURE 8 shows the creep displacement at the studied temperature.

FIGURE 7 -Time to failure as function of the temperature for creep tests in tension at 50% of T ult FIGURE 8 -Normal displacement depending on time during creep tests at 50% of T ult at 20°C, 35°C and 50°C

Modelling approach

To describe the observed creep behaviour, the Burger model is used [START_REF] Majda | A modified creep model of epoxy adhesive at ambient temperature[END_REF]. This model gives the creep strain as a function of the time , the creep stress and three parameters , and with the following equation:

The creep parameters were identified for all the presented Arcan tests. It enables to develop a nonlinear Burger model for which the parameters depend on the creep stress.

To verify the accuracy of this modelling approach, the nonlinear Burger law was implemented in a numerical model of the fastener C-Claw in the case of tension load. The used model is the same than the one in the study of Djeumen et al. [START_REF] Djeumen | Creep investigations on adhesively bonded fasteners developed for offshore steel structures[END_REF] and Loiseau et al. [START_REF] Loiseau | Inverse identification method of adhesive creep properties from real scale investigations on bonded fastener[END_REF]. The results of these simulations were then compared to the experimental data from the investigations at the scale of the fastener presented in the study of Loiseau et al. [START_REF] Loiseau | Determination of Creep Behaviour of Adhesively Bonded Assembly -Application to Adhesively Bonded Steel Fasteners[END_REF]. FIGURE 9 compares the results of the numerical simulation with the experimental data for the creep investigations at 60% of the ultimate capacity in tension of the fastener. We can notice a good correlation between the numerical model and the experimental data.

FIGURE 9 -Comparison of the numerical simulation with the experimental data of the creep tensile tests at 60% of the ultimate capacity in tension of the fastener

Conclusion

In this study, experimental investigations of the creep behaviour have been carried out at the scale of adhesive joints using Arcan test. To do so, dedicated resources were developed: a creep mechanical frame and a specific measurement system for the local displacements.

The Burger model was used to describe the observed creep behaviour. In tension, three load levels and three temperatures were investigated. This enabled to develop a nonlinear creep law. Consequently, the dependency of the creep parameters with the stress is considered.

The nonlinear model was implemented in a numerical model (finite element) of the fastener. The results were compared to experimental data of creep tests performed at the scale of the fastener in previous studies. A good correlation is observed between the numerical and the experimental data.

This offers interesting perspectives to develop a creep model from experimental investigations at the scale of the adhesive joint to predict the creep behaviour of bonded assemblies at the scale of the structure.

FIGURE 2 -

 2 FIGURE 2 -(a)Arcan test TCS developed by Créac'hcadec et al.[START_REF] Créac'hcadec | A stress concentration-free bonded arcan tensile compression shear test specimen for the evaluation of adhesive mechanical response[END_REF], (b) Zoom on the Arcan TCS beaks[START_REF] Créac'hcadec | A stress concentration-free bonded arcan tensile compression shear test specimen for the evaluation of adhesive mechanical response[END_REF]