Thomas Miro Padovani

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

Ahmed Ketfi-Cherif
email: ahmed.ketfi-cherif@renault.com

Yann Chamaillard
email: yann.chamaillard@univ-orleans.fr

Yann Chamaillard Optimal Energy

Miro - Padovani

Optimal Energy Management Strategy for a Hybrid Vehicle provided with a Dual Electric Storage System

Keywords: Energy management strategy, Optimal control, Hybrid electric vehicle, Dual storage system, Pontryagin's Minimum Principle

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

Optimal Energy Management Strategy for a Hybrid Vehicle provided with a Dual Electric Storage System

Thomas MIRO PADOVANI * Guillaume COLIN * * Ahmed KETFI-CHERIF * Yann CHAMAILLARD * *

INTRODUCTION

Hybrid electric vehicles (HEVs) and Plug-in hybrid electric vehicles (PHEVs) are currently among the most promising solutions to maximize global powertrain efficiency and therefore reduce the fuel consumption of our vehicles. Powertrain hybridization opens the way to an almost unlimited set of solutions in terms of architectures, technologies and sizing, each leading to a different compromise. The main feature of the powertrain discussed in this paper comes from its electric storage system. While the electric storage system of most HEVs is made of a single Li-ion or Ni-MH battery pack, this paper focuses on a singular architecture, composed of a 12V Li-ion battery associated in series with a pack of double layer capacitors (DLC). Energy transfer between the two devices is performed by a bi-directional DC/DC converter. The main motivation behind this architecture is to combine the high specific power and durability of the DLC with the high specific energy of Li-ion batteries. With an appropriate supervisory control, this system can lead to reduced energy losses while enhancing battery life span. The work in this paper aims at proposing an embeddable energy management strategy (EMS) for the studied Mild-HEV, based on optimal control theory, and taking into account both storage devices.

Optimal control theory has been widely investigated to provide both off-line and on-line oriented solutions for the HEV energy management problem. Off-line strategies such as deterministic dynamic programming (DDP) [START_REF] Sunström | Optimal hybridization in two parallel hybrid electric vehicles using dynamic programming[END_REF]; [START_REF] Lin | Power management strategy for a parallel hybrid electric truck[END_REF]; [START_REF] Elbert | Implementation of dynamic programming for n-dimensional optimal control problems with final state constraints[END_REF] and Pontryagin's Minimum Principle (PMP) require the knowledge of future driving conditions to achieve optimality.

On the other hand, on-line strategies need to be causal, and while they can be implemented in a vehicle, optimality is no longer guaranteed. The most common optimalcontrol based on-line strategies are stochastic dynamic programming (SDP) [START_REF] Moura | A stochastic optimal control approach for power management in plug-in hybrid electric vehicles[END_REF]; [START_REF] Kolmanovsky | Optimization of powertrain operating policy for feasibility assessment and calibration: Stochastic dynamic programming approach[END_REF] and the Equivalent Consumption Minimization Strategy (ECMS) [START_REF] Paganelli | Equivalent consumption minimization strategy for parallel hybrid powertrains[END_REF]; [START_REF] Chasse | Online optimal control of a parallel hybrid with costate adaptation rule[END_REF], which is the on-line counterpart of the PMP.

Dual storage systems embedded in vehicular applications have already been addressed in previous work [START_REF] Romaus | Optimal energy management for a hybrid energy storage system combining batteries and double layer capacitors[END_REF]; [START_REF] Allegre | Infuence of control strategies on battery/supercapacitor hybrid energy storage systems for traction applications[END_REF]; [START_REF] Cheng | Active control of power sharing in a battery/ultracapacitor hybrid source[END_REF]; [START_REF] Omar | Power and life enhancement of battery-electrical double layer capacitor for hybrid electric and charge-depleting plug-in vehicle applications[END_REF]. The work presented in this paper differs from the singular electric architecture considered where the 12V battery is a part of both the 12V and the high voltage (HV) traction grid. Moreover, the optimization of the dual storage system is encompassed in the traditional HEV energy management problem by considering an extra state variable similarly to [START_REF] Fontaine | Optimal control of a parallel hybrid vehicle equipped with a dual electrical storage system[END_REF]. Both the off-line optimal solution and an on-line oriented solution are presented and discussed. The powertrain architecture and models are presented in section 2. The optimal control problem dealing with both storage devices is presented in section 3. In section 4, optimal results yielded through DDP are discussed. Finally, an application of the ECMS for the dual storage system is proposed in section 5, and the results are discussed and compared to the DDP optimal solution.

INVERTER

DC/ DC

Battery

VEHICLE SPECIFICATIONS AND MODELS

Vehicle architecture

The vehicle considered in this paper is a parallel Mild-HEV, provided with a 12 kW belt-driven starter generator connected to the crankshaft of a 1.5 L turbo diesel engine producing 80 kW. The engine is connected to a 6 gear manual transmission. The main innovation in the powertrain comes from its electric storage system, composed of two storage devices connected in series: a 12V Li-ion battery pack and a 36V DLC pack. Both storage devices in series make the high voltage (HV) grid (48V), therefore both devices contribute to traction power. The 12V battery pack is connected to the HV grid via a bi-directional DC/DC converter. Finally, the inverter supplying the electric machine (EM) is connected to the HV grid. Both storage devices offer a combined usable electric energy of approximately 200 W.h. Fig. 1 illustrates the electrical diagram of the architecture as well as the naming convention for currents and voltages in the system. The singularity of this topology, where the 12V battery belongs to both the 12V grid and the HV grid, logically results in a strong coupling between the two grids of the vehicle. Consequently, the power consumption of the inverter P IN V will directly impact the 12V grid through the battery; similarly the power consumption P A of the auxiliaries on the 12V grid will impact the available power for the EM. The main powertrain and vehicle characteristics are listed in table 1.

Powertrain and vehicle models

The model used for simulation is a control-oriented backward quasi-static model. The speed profile is given by a driving cycle which provides both the vehicle speed v(t) and acceleration setpoint. The vehicle model computes the corresponding required powertrain torque at the wheel T sp pwt as follows:

T sp pwt (t) = (m veh . v(t) + F res (t) + F brk (t) +m veh .g. sin(α(t))).r wheel , (1)
with m veh the vehicle mass, F brk the force applied by the braking system, α the slope angle, g the gravitational acceleration, and r wheel the wheel radius. F res is the resistive force applied to the vehicle caused by the rolling friction and drag. It is modeled by a second order polynomial:

F res (t) = p 1 .v 2 (t) + p 2 .v(t) + p 3 , (2)
where p 1-3 are identified via the experimental road load given by the deceleration law of the considered vehicle.

The manual transmission model provides the rotary speed of the internal combustion engine (ICE) and EM, ω ICE and ω EM , given the wheel speed ω wheel :

ω ICE = Rgb(g nbr ).ω wheel ω EM = R EM .ω ICE , (3)
where Rgb(g nbr ) expresses the gearbox ratio corresponding to the gear number g nbr and R EM is the constant gear ratio between the crankshaft and EM.

The engine consumption ṁfuel is computed using a lookup table given by experimental results on an engine test bench in nominal operating conditions:

ṁfuel = f ICE (ω ICE , T ICE ). ( 4
)
The dynamics on the torque are ignored, and T ICE the output torque of the ICE exactly matches the torque setpoint given by the EMS.

The EM electrical power P IN V is given by a look-up table based on dynamometer testing,

P IN V = f EM (ω EM , T EM ), (5)
where T EM is the output torque of the EM. As for the ICE, the EM's dynamics is ignored and T EM matches a given valid setpoint.

To ensure that the driver torque demand is always met, the EM torque and ICE torque are linked by the following:

T pwt = (T ICE + T EM .R EM ).Rgb(g nbr ) = T sp pwt . (6)
The battery and DLC pack are modeled as an open circuit voltage (OCV) source associated with an internal resistance R in series, both depending on the battery state of energy (SoE):

S v : V 1 = OCV 1 (SoE) + R 1 (SoE).I 1 V 2 = OCV 2 (SoE) + R 2 (SoE).I 2 V HV = V 1 + V 2 , (7)
with V 1 and V 2 the battery and DLC voltages. As a convention, the index 1 refers to battery-related quantities, and 2 refers to the DLC. The bi-directional DC/DC converter is modeled as follows:

P C = V 1 .I 4 = η C .V HV .I 3 , (8)
with

η C = η if P C ≥ 0 1 η else, ( 9
)
where η is a constant approximation of the DC/DC efficiency. The auxiliaries power consumption P A is considered as a known constant perturbation. Finally, the electrified part of the powertrain presented in Fig. 1 can be modeled from Kirchhoff's laws with the following system of equations:

S :          V 1 .I 4 = η C .V HV .I 3 ; I A = P A V 1 I IN V = P IN V V HV ; P C = V 1 .I 4 I 2 = I IN V -I 3 ; I 1 = I 4 + I 2 -I A , (10)
where as a convention, a positive power or current will result in charging the storage devices. S can be reduced to a system with 3 equations and 3 unknowns, namely I 3 , I 2 and I 1 :

S I : I 3 = P C η C .V HV ; I 2 = P IN V V HV -I 3 ; I 1 = P C -P A V 1 + I 2 .
(11) When incorporating the expressions of the voltages in S v into the system of equations S I , one obtains a nonlinear system of second degree equations with no simple analytical solution. The difficulty of the system comes from the coupling between voltages and currents, coming from the internal resistance of the storage devices. In order to break the algebraic constraint between the computation of currents and voltages, a fixed point method is proposed. For the first step, we assume as initial conditions R 1 = R 2 = 0, allowing us to replace V 1 and V HV in S I by respectively OCV 1 and OCV 1 + OCV 2 . This assumption leads to a first estimation of the currents in S I . These current estimations are used to compute the voltages in S v , which are in turn used to compute a better estimation of the currents in S I . By repeating the process, one iteratively computes a more precise estimation of the currents in S I . Once the currents converge with a successive error below 1%, the corresponding solution is selected. For the electrical system considered here, a maximum of 5 iterations was enough to meet the convergence criterion of a 1% error or less. Once the 3 currents of S I are known, it is possible to calculate the dynamics of both storage devices' SoE by integrating their inner power P inn over time:

ṠoE j (SoE j , t) = P inn j E j = OCV j (SoE).I j (t) E j , (12)
where j = {1, 2} depending on the considered device, E j is the battery or DLC available energy when fully charged, and P inn j is the inner battery or DLC power.

OPTIMAL CONTROL PROBLEM FORMULATION

Optimal control based EMS require the definition of an optimality criterion, also called performance index, which serves as the control policy's objective. The performance index considered in this paper is the global fuel consumption along the trip:

J(x(t), u(t)) = t f t0 ṁfuel (u(t), t).dt + Φ (x(t f )) u opt = argmin u J(x(t), u(t)) (13) with Φ (x(t f )) = 0 if x(t f ) ≥ SoE target ∞ else. ( 14
)
The control variable u is chosen as follows:

u = u 1 u 2 = P C P IN V . ( 15
)
The considered state variable x is made of both storage devices' SoE:

x = x 1 x 2 = SoE 1 SoE 2 , (16)
Φ (x(t f )) is a function ensuring that the solution meets the final requirement on both SoEs. The final SoEs, namely SoE target , will be chosen equal to their initial values to ensure sustaining operation. The optimization problem ( 13) is subject to :

     u(t) ∈ U (t) x(t) ∈ χ ẋ1 (t) = f 1 (x 1 (t), u(t), t) , x 1 (t 0 ) = x 1,0 ẋ2 (t) = f 2 (x 2 (t), u(t), t) , x 2 (t 0 ) = x 2,0 , (17)
where U (t) is the set of admissible commands for P IN V and P C . Likewise, χ is the space of admissible values for both state variables.

The two state functions f 1 and f 2 are given by :

f 1 (x 1 (t), u(t), t) = OCV 1 (x 1 (t)).I 1 (x 1 (t), u(t), t) E 1 , f 2 (x 2 (t), u(t), t) = OCV 2 (x 2 (t)).I 2 (x 2 (t), u(t), t) E 2 . ( 18
)
with the expression of I 1 and I 2 given by ( 11).

One of the distinctive features of the optimal control problem presented lies in the calculation of the set of admissible commands U (t). In addition to regular limitations on the actuators:

P IN V ≤ P IN V ≤ P IN V P C ≤ P C ≤ P C , (19)
where X and X denote respectively the minimum and maximum admissible values for X. Battery and DLC power limitations must be taken into consideration:

P 1 ≤ P 1 ≤ P 1 P 2 ≤ P 2 ≤ P 2 , (20)
Battery and DLC power limitations are needed to ensure that the optimal commands yielded by the EMS will be accepted by the battery management system (BMS), and therefore will not deteriorate the storage devices. Power limitations may also be used to ensure a set of voltage limitations, for instance on the 12V grid. As a result, 4 inequality constraints are considered to determine the final set of admissible commands U (t) :

             P C ≤ P C ≤ P C P IN V ≤ P IN V ≤ P IN V P 2 .C 2 ≤ P IN V - 1 η C .P C ≤ P 2 .C 2 (P 1 + P A ).C 1 ≤ P IN V + (C 1 - 1 η C ).P C ≤ (P 1 + P A ).C 1 , (21) with
C 1 = V HV V 1 , C 2 = V HV V 2 . ( 22
)
The inequality constraints corresponding respectively to the DC/DC, inverter, DLC, and battery power limitations are expressed via ( 11) and ( 19). The inequalities system (21) takes the general formulation:

     u 1 ≤ u 1 ≤ u 1 u 2 ≤ u 2 ≤ u 2 z ≤ u 2 + α.u 1 ≤ z w ≤ u 2 + β.u 1 ≤ w, (23)
with in this case α < 0 and β > 0. From (23), one can determine global limitations on u 2 versus u 1 : max(u 2 , z-α.u 1 , w-β.u 1 ) ≤ u 2 ≤ min(u 2 , z-α.u 1 , w-β.u 1 ) (24) Eq. ( 24) gives a set of 9 inequalities, 6 of which are relevant to define the absolute limitations on u 1 :

max(u 1 , w -u 2 β , z -u 2 α , w -z β -α ) ≤ u 1 ≤ min(u 1 , w -u 2 β , z -u 2 α , z -w α -β ).
(25) Finally, with (25) giving the global limitations on u 1 ensuring the set of constraints (23), and with (24) giving the global limitations on u 2 versus u 1 , the set of admissible commands U (t) ∈ (u 1 , u 2 ) is defined. Fig. 2 plots U (t) for nominal operating conditions. The numerical values considered for this example are the following:

-10000 ≤ P IN V ≤ 10000, -1000 ≤ P C ≤ 2000 -2500 ≤ P 1 ≤ 2500, -8000 ≤ P 2 ≤ 8000 (26)

The domain corresponding to the actuators' limitations is represented by the blue lines. The green dashed lines represent the controls reaching the minimum and maximum battery power limitations, whereas DLC power limitations are represented in black. The set of admissible commands is the intersection of all previously mentioned domains, and is represented in red. Determining a priori the space of admissible controls rather than eliminating non-feasible commands a posteriori has two main advantages. Not only does it ensure feasible commands, but it also reduces the search area of the control space to what is necessary only during the minimization process, reducing the computational burden for vehicle embedded applications.

OFF-LINE OPTIMAL SOLUTION

The optimal control problem (13) presented in section 3 is solved using DDP. Commonly, the EMS problem only considers one state variable corresponding to the single storage device present. The extra state variable considered in this paper greatly increases the computational effort required by the DDP to solve the problem. However, computational time remain reasonable for off-line computer based simulations; mainly because the very low levels of battery energy considered require a smaller number of points to properly discretize the state space. Fig. 3 and 4 present the optimal solution for the Worldwide harmonized Light duty driving Test Cycle (WLTC) considering 250 W of auxiliaries' power consumption and charge sustaining operation. When analyzing the optimal control policy, several interesting behaviors can be identified. First of all regarding the DC/DC converter, the reversibility of the latter is never used, as swapping energy back and forth between storage devices leads to unnecessary losses when all the driving conditions are perfectly known a priori. Therefore the DC/DC converter is only used to charge the 12V battery, strictly to compensate for the auxiliaries power consumption. The optimal control policy of the DC/DC converter shows that the instantaneous power P C (t) is either null or close to the maximum efficiency of the converter. However, the mean value of P C (t)) along the entire driving cycle is always equal to the mean value of the auxiliaries' consumption.

As for the P IN V command, it is worth noting that for this specific powertrain, using the EM to optimize the ICE's operating point is not optimal. The EM power and efficiency are not high enough to justify operating point optimization for the diesel engine which has a large area of flat specific fuel consumption. As a result, the EM is used as a generator only when the electric energy recovered during deceleration phases is not enough to compensate for the auxiliaries' consumption.

When analyzing the impact of the optimal control policy on the storage devices, it appears clearly that the EMS avoids high battery currents. Indeed, during deceleration phases when the EM produces its maximum levels of power, the DC/DC power goes back to 0 to avoid overcharging the battery, cf Finally, when analyzing both SoEs, the DLC clearly serve as the main energy buffer, whereas SoE 1 remains almost constant, as shown in Fig. 4. Moreover, SoE 2 is highly correlated to vehicle speed: the optimal EMS always ensures depletion between two deceleration phases to guarantee sufficient available energy storage for each deceleration phase. Therefore, SoE 2 is in average higher at low speed than it is at high speed.

PROPOSED ON-LINE STRATEGY

To solve the optimal control problem (13) on-line, an extended version of the well-known ECMS taking into account two state variables is proposed. The PMP states that the optimal policy u opt (t) has to minimize, for all t ∈ [t 0 ; t f ], the following Hamiltonian function:

H(u, λ, x, t) = ṁfuel (u 1 , t) + λ 1 (t). ẋ1 (u, x 1 , t) +λ 2 (t). ẋ2 (u, x 2 , t) H eq (u, s, x, t) = ṁfuel (u 1 , t).H LHV -s 1 (t).P inn 1 (u, x 1 , t) -s 2 (t).P inn 2 (u, x 2 , t), (27)
with

s(t) = s 1 (t) s 2 (t) =      -λ 1 (t).H LHV E 1 -λ 2 (t).H LHV E 2      , (28)
where .) ∂x is close to 0 for a large range of SoE. This assumption is verified here by calculating the optimal values of s 1 (t) and s 2 (t) using DDP, cf Fig. 5, as the costate of the dual problem is given by:

λ(t) = λ 1 (t) λ 2 (t
λ(t) = ∂ ∂x V opt (x(t), t), (30)
where V opt (.) is the optimal-value, or cost-to-go function considered in Bellman's principle [START_REF] Elbert | Implementation of dynamic programming for n-dimensional optimal control problems with final state constraints[END_REF].

The simulation results presented in this section were obtained considering an on-line adaptation of s 1 (t) and s 2 (t) via a proportional-integral controller of the SoE for each:

s j (t) = s j 0 + K j p .SoE j err + K j i . where j = {1, 2} and the parameters s j 0 , K j p , and K j i are kept unique for every driving cycle. Therefore they are chosen in order to obtain a global robustness on the SoE regulation for all driving conditions, especially to ensure that the SoEs remain bounded within the state constraints domain χ. Readers can refer to [START_REF] Chasse | Online optimal control of a parallel hybrid with costate adaptation rule[END_REF] for more information on equivalence factor tuning methods. While an on-line adaptation of s 1 and s 2 using a PI controller is by nature sub-optimal, it is essential to prevent both SoEs from diverging when driving conditions are not known a priori. Fig. 6 presents simulation results given by the proposed on-line strategy on the WLTC cycle with P A = 250W. Several differences appear regarding the optimal solution presented Fig 3 and4. While the DC/DC is mainly used to charge the 12V battery, here it is also used in reverse mode a few times to ensure robustness of the SoE 2 regulation. Moreover, SoE 2 is no longer correlated with vehicle speed, as the proposed on-line strategy does not anticipate future acceleration or deceleration phases. Apart from these items, all of the behaviors discussed in section 4 are common to the optimal solution and the on-line strategy presented.

The sub-optimality results of the on-line strategy are given in table 2 for different driving cycles and accessories consumption. The sub-optimality is defined as the relative extra fuel consumption regarding the optimal DDP solu- tion. The driving cycles considered are the WLTC and the Artemis cycles. The cycles were repeated ten times in a row, in order to ensure that the possible energy unbalance between x(t 0 ) and x(t f ) of a few W.h was absolutely negligible regarding the total fuel mass consumed. Overall, the extra fuel consumption yielded by the on-line strategy did not exceed 2.8%, it even remained below 1% for P A = 250 W. But most importantly, these results were obtained with a single set of tuning parameters for the SoE regulators, highlighting the potential of the ECMS strategy to give near-optimal results with a good overall robustness and a low tuning effort.

CONCLUSION

Dual storage systems are viable solutions for electrified powertrains to combine the advantages of different technologies or to benefit from the energy stored in the traditional 12V battery. The system presented in this paper combines a 12V Li-ion battery with a DLC pack, both devices in series forming the HV traction grid. Energy transfer between the 12V and HV grids is ensured by a bidirectional DC/DC converter. The electric system is modeled and the algebraic constraint resulting in the coupling of voltages and currents is solved through a fixed point method. The energy management of the HEV powertrain is formulated as an optimal control problem considering two control variables and two state variables. A mathematical approach is presented to define the space of admissible controls taking into account actuators and storage devices power limitations. The optimal control problem is solved off-line using DDP, and the main results of the optimal solution are discussed. Although the optimal solution's only target is to minimize fuel consumption, it also preserves the battery SoH by avoiding high currents. High currents are the main cause of energy losses via the joule effect, and at the same time a factor aggravating battery aging.

Finally, an on-line EMS based on the ECMS is presented, taking into account the two state variables with a double co-state formulation. The strategy shows very promising simulation results, yielding a fuel consumption very close to the optimal solution for the considered system in all experimented situations. The overall robustness of the strategy is one of its main assets, as all the simulations were conducted with a single set of tuning parameters.

Fig. 2 .

 2 Fig.2. Depiction of the 2 dimensional space of admissible commands U for the considered optimal control problem.

  Fig 3. Symmetrically, during engine starts where the EM requires a lot of power, the DC/DC converter is used to relieve the battery by lowering its discharge current. The EMS always favors using the DC/DC to avoid too high values for |I 1 |, even if it means increasing |I 2 |. This comes from the fact that R 1 is substantially higher than R 2 , the battery is therefore more subject to ohmic losses. Avoiding high battery current is also very important to preserve the battery state of Health (SoH), as it also reduces battery temperature and overall A.h throughput, which are the three main factors aggravating battery aging Miro-Padovani et al. (2013).

Fig. 3 .Fig. 4 .

 34 Fig.3. Zoom on optimal control trajectory during deceleration phases. When P IN V is high, P C remains null to avoid increasing the battery current.

Table 1 .

 1 Vehicle characteristics.

	𝑰 𝑰𝑵𝑽	𝑰 𝟐
	𝑰 𝟑	DLC		𝑽 𝟐
	(HV grid) 𝑽 𝑯𝑽			𝑰 𝑨
		𝑰 𝟒	𝑰 𝟏
				A
				𝑽 𝟏 (12V grid)	U X
		LTO
	Fig. 1. Electrical diagram of the dual storage system
	considered.
	Vehicle mass		1500 kg
	Engine displacement		1.5 L
	Engine maximum power		80 kW
	Electric machine maximum power	12 kW
	Combined energy storage capacity 200 Wh

Table 2 .

 2 Fig. 6. On-line strategy results on the WLTC cycle. Extra fuel consumption of the on-line strategy relative to the optimal solution.

	Speed [km/h]	50 100
				0
				90						SOE 2
	SoE [%]	60 70 80						SOE 1
				50
				2.6 2.8						s	2
	s [-]	2.2 2.4						s 1
				2
				1.8
	P C [W]	0 1000
	-1000
	10000
	P INV [W]	0 5000
	-5000
				0	200	400	600	800	1000	1200	1400	1600
								Time [s]
					Extra fuel consumption [%]
	P A [W]	WLTC	Urban Extra-urban Traffic jam
	250		0.7		0.4		0.5	0.4
	600		1		1.3		0.7	1.7
	900		1.2		1.8		0.8	2.8


Antoine Bérut

Artyom Petrosyan

Juan Ruben Gomez-Solano

Sergio Ciliberto

Fluctuations in an aging system: the

Keywords: Function

Introduction and Motivations 1.Gelatin and the sol-gel transition

Gelatin is a thermoreversible gel [START_REF] Nijenhuis | Thermoreversible networks: viscoelastic properties and structure of gels[END_REF]. It is a heterogeneous mixture of water-soluble denatured collagen protein chains, extracted by boiling animal by-products (skin, tendons, ligaments, bones, etc.) in water. Collagen molecules are rods of 300 nm length, made of three strands, with high average molecular weights. This triple-helix structure is stabilised by hydrogen bonds and has a diameter of ∼ 1.4 nm. The chemical treatment used to produce gelatin breaks crosslinks between strands, but can also hydrolyze strands into fragments. Thus a broad molecular weight distribution is obtained for gelatin [START_REF] Djabourov | Gelation of aqueous gelatin solutions. i. structural investigation[END_REF][START_REF] Djabourov | Gelation of aqueous gelatin solutions. ii. rheology of the sol-gel transition[END_REF].

Above a temperature T melt ∼ 40 • C, the gelatin chains are in coil conformation. The gelatin solution is in a viscous liquid phase, called "sol" phase. Below a temperature T gel ∼ 30 • C renaturation of the native triple helix structure occurs, and chains form a percolating three-dimensional network of helical segments connected by single strand coils. The gelatin solution is in an arrested state with elastic behaviur, called "gel" phase. The coil-helix transition is completely reversible and the transition from one phase to the other is called the "sol-gel" transition [START_REF] Nijenhuis | Thermoreversible networks: viscoelastic properties and structure of gels[END_REF].

Physical properties of the sol phase, and of the sol-gel transition are studied in [START_REF] Bohidar | Kinetics of sol-gel transition in thermoreversible gelation of gelatin[END_REF][START_REF] Bohidar | Study of sol-state properties of aqueous gelatin solutions[END_REF] for different gelatin concentrations above 4 wt%. In particular, it was seen that there are at least three successive steps in the transition: monomer to aggregate formation, random-coil-single-helix transition (disorder-order transition), and single-helix-triple-helix transition (order-order transition). It is then possible to identify different phase states in the sol domain: the sol state I where the chains have random coil conformations and the sol state II where single and triple helices begin to form (without reaching gelation).

The gel phase was also shown to share properties with glassy materials, which are outof-equilibrium metastable systems. After a quench at T < T gel the system is frustrated by topological constraints because each gelatin chain is involved in at least two helices, and neighbouring helices are competing for the shared portions of non-helical chain. Therefore, the system displays physical ageing: its physical properties slowly evolve with time, through a process known as structural recovery. For example, the small-strain shear modulus of a 5 wt% gelatin solution quenched at 20 • C increases logarithmically as a function of the ageing time [START_REF] Ronsin | Interplay between shear loading and structural aging in a physical gelatin gel[END_REF]. And the elasticity of gelatin gels during slow cool and heat cycles exhibits memory and rejuvenation effects similar to the ones found in spin glasses [START_REF] Parker | Glassy dynamics of gelatin gels[END_REF].

Although it is known that mechanical properties of gelatin gels are very sensitive to temperature variations, previous thermal history of the gel, and time, this system has some interesting experimental features:

• The fact that the transition is thermoreversible allows us to do melting/gelation cycles simply by controlling the temperature of the sample.

• The ageing rate can in theory be controlled by changing the quench depth.

• The length-scale of the collagen chains (300 nm) is big enough to be sensed by a micro-particle of 2 µm.

This particular sol-gel transition was chosen for previous works done in our group about fluctuations of Brownian particles in quenched gelatin samples [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF].

A summary of the previous works results and our motivations are presented in the next section.

Previous work: anomalous variance, heat flux and Fluctuation Dissipation Theorem violation in an ageing bath

Previous works [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF] showed that a particle trapped with optical tweezers in a liquid droplet of gelatin solution, quenched at a temperature below T gel exhibits anomalously high position fluctuations right after the quench. These anomalous fluctuations could be interpreted as an effective temperature, and were consistent with a violation of the Fluctuation Dissipation Theorem (FDT) and with an exchange Fluctuation Theorem (xFT) for the heat exchanged between two heat baths at different temperatures [START_REF] Jarzynski | Classical and quantum fluctuation theorems for heat exchange[END_REF].

We reproduce here some figures from [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF] and recall the associated key results:

• The variance of the position σ 2 x = x 2 exhibits anomalously high value for short times (∼ 5 s) right after the quench. It then stabilises at the equipartition value σ 2

x eq = k B T /k for ∼ 200 s (where k B is the Boltzmann constant and k the stiffness of the trap). And it finally decreases logarithmically for long times after the quench. See figure 1a. The anomalously high variances can be interpreted in terms of effective temperatures:

σ 2 x = k B T eff /k > σ 2
x eq . • The Probability Distribution Functions of position fluctuations are Gaussian at any time after the quench, but their variances decreases with time (in agreement with the variances observed). See figure 1b.

• The Probability Distribution Function of the heat exchanged between the particle and the bath during short times after the quench is asymmetrical. See figure 2.

• The asymmetry function of the heat ρ(q) = ln (P (q)/P (-q)) (where P (q) is the probability of observing the value q of the heat) satisfies an exchange Fluctuation Theorem: ρ(q) = ∆βq, with its slope ∆β directly linked to the effective temperatures defined from the variances (see [START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF]).

• The Fluctuation Dissipation Theorem is violated only for short times after the quench, and this violation can be linked with the amount of heat exchanged between the particle and the bath during the same time. See figure 3.

Unfortunately, none of those key results was found to be reproducible, and we believe that they were only due to an artifact in the data and/or in the analysis method. Therefore, we present in this article a detailed and careful analysis of trajectories of particles trapped in a droplet of gelatin solution quenched at a temperature below T gel . We show that there is indeed no effective temperature for this system, which surprisingly exhibits equilibrium-like properties while aging. The article is organized as follow: in the first section we describe our experimental setup, the gelatin solution preparation, optical trapping system, and local heating for fast quenching method. In the second section, we discuss our experimental results: we show with some bulk measurements that the system is indeed aging after the quench, we then analyze the effect of a slow drift in trajectory and of a mixing between time and ensemble averages in data analysis to show that there is no anomalous variance of the particle's position. We end by discussing the consequences of the absence of effective temperature for the heat exchange, Fluctuation Theorem and Fluctuation Dissipation violation.

2 Experimental set-up

Gelatin sample preparation

We use gelatin powder from porcine skin, produced by Sigma-Aldrich : gel strength ∼ 300 g Bloom, Type A, BioReagent, suitable for cell culture. This gelatin is derived from acid-cured tissue, whereas type B is derived from lime-cured tissue. We work with gelatin at a weight concentration of 5 wt%. The samples are prepared following a standard protocol [START_REF] Normand | Gelation kinetics of gelatin: A master curve and network modeling[END_REF]: the wanted amount of powder is dissolved in bidistilled water, which is then heated for ∼ 30 min at ∼ 60 • C while slowly stirred until the solution is transparent and homogeneous. While the solution is still liquid, ∼ 2 mL are filtered using a Millex syringe driven filter unit with 0.45 µm pore size mixed cellulose esters membrane. Then 15 µL of an aqueous solution of silica beads (radius R = 1.00 ± 0.05 µm) with concentration 10 7 particle mL -1 are added, and the solution is strongly agitated. The non-filtered and final solutions are let gel at room temperature and kept in the refrigerator for later use.

We use a disk-shaped glass cell, with an Indium Tin Oxyde (ITO) coated microscope slide, a free-volume to avoid problems if the volume of solution changes during gelation, and a Wavelength Electronics TCS10K5 thermal sensor for temperature measurement (see figure 4). To fill the cell, the gelatin solution with dispersed silica particles is taken from the refrigerator and heated at ∼ 50 • C until it is in the sol phase.

Optical trapping and controlled gelation

For different purposes we used two different alignments of an optical tweezers set-up:

The first set-up uses a laser beam (λ = 532 nm) separated in two cross-polarized beams which enables us to have two traps with no interference between them. A laser diode (λ = 980 nm) is aligned with the green laser and used to heat locally the sample. The tracking is done using a fast camera which is able to track two particles at 600 Hz. The acquisition rate is not very fast but by previously using a calibration target, we can directly convert the displacement of the particle from pixels to µm. See figure 5.

The second set-up uses a single laser diode (λ = 980 nm) to trap the particle and to heat locally the sample. A He-Ne laser (λ = 632.8 nm) is aligned with the laser diode and deflected by the particle. This deflection is measured using a position sensing diode which is able to track one particle at more than 10 kHz. Contrary to the first set-up, the position signal is in arbitrary unit and a supplementary calibration is needed for each measurement to convert the trajectory of the particle in physical units. See figure 6.

For both set-ups, the microscope objective is an oil-immersion Leica HCX PL. APO ×63 with high numerical aperture N.A. = 1.4. The microscope objective is surrounded with a custom-made heating ring, made with a Minco flexible resistor and a Wavelength Electronics TCS10K5 thermal sensor for temperature measurement. A feedback control is managed by the temperature control module (TCM-39032) of a modular laser diode controller (ILX Lightwave LDC-3900). As already mentioned, another thermal sensor is inserted directly inside the cell (see figure 4) and another feedback control is done by mastering the current going through the ITO-coated microscope slide, with an Instec MK1 Board and a PID software. These two temperature devices ensure that the temperature of both the microscope objective and the cell are well controlled. The precision achieved on the temperature control is about ±0.05 • C.

To trap particles, the temperature of both the microscope objective and the cell are set to 38 • C so that the gelatin is in the sol phase. Then, one particle is found and trapped at a given distance from the bottom surface of the cell (typically h = 15 µm). Then, the temperature controls are set to a value below the gel transition (typically T fb = 27 • C) and we let the sample gel for a few hours (typically between 6 and 10 hours). Since the gelatin shows a lot of hysteretic behaviour [START_REF] Djabourov | Gelation of aqueous gelatin solutions. i. structural investigation[END_REF][START_REF] Djabourov | Gelation of aqueous gelatin solutions. ii. rheology of the sol-gel transition[END_REF], this gelation procedure appears to be important and the rheology of the gel can vary if the gelling time is very different (e.g. a few days). Moreover, one must pay attention to regularly check the distance between the bead and the bottom surface of the cell, since the focal distance of the microscope objective always drifts slowly when its temperature is changed.

The refractive index of the liquid gelatin solution was measured n gel = 1.3415, which An infrared laser diode is used to trap and to melt the gelatin. The deflection of a He-Ne laser induced by the trapped bead is measured by the position sensing diode (PSD). The white light source and camera are only used for direct visualization but not for the measurements. "BS" is a beam-splitter, "IF" and interferential filter to suppress the infrared beam and "DM" are dichroic mirrors. is close to value in water n water = 1.3335. It follows that the trapping stiffnesses in gelatin solution should be close to the ones in water with same experimental parameters.

Local heating and fast quenching method

When the particle is trapped and the sample is properly gelled at a given controlled temperature T fb < T gel , the local quenches are done in a similar way than the one presented in [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF]: The power of the 980 nm laser diode is risen to a high value (typically1 230 mW) during a given time (typically τ melt = 200 s). Because of the light absorption of the water molecules in the solution, the temperature of the gelatin around the particle (which is at the focal point of the microscope objective) increases by a small amount δT . Following the formula in reference [START_REF] Peterman | Laser-induced heating in optical traps[END_REF]:

δT = P α 2πK ln 2πh λ -1 (1)
where α = 50 m -1 is the attenuation coefficient of water at 27 • C for wavelength 980 nm [START_REF] Palmer | Optical properties of water in the near infrared[END_REF], and K = 0.61 Wm -1 K -1 is the thermal conductivity of water. Here, we await:

δT 11 • C. (2)
This increase in temperature is only roughly estimated. Especially because we don't really know what is the absorption of the microscope objective for the near infrared, and because it is impossible to measure the temperature with a usual probe on this very small scale. But it is seen that the increase is strong enough to melt a small droplet of gelatin (radius R d ∼ 10 µm) around the bead. Then, the power is quickly decreased to a low value (in the case where the same laser diode is used to trap and heat) or to zero (in the case where another laser is used to trap the particle), and the sample is let gel for a given time (typically τ rest = 500 s). Since the thermal diffusivity of water is κ = 0.143 × 10 -6 m2 s -1 at 25 • C, the time τ κ needed to dissipate the heat from the droplet to the bulk is short:

τ κ ∼ R 2 d κ ∼ 2 × 10 -4 s. (3)
Hence, the gelatin is believed to experience a fast quench at temperature T fb < T gel and should start ageing 2 . After the resting time τ rest at low temperature T fb , the power of the laser diode is risen again, and another quench is done. Note that the exact duration of τ rest was not considered as important, because it was believed that the melting "resets" the gelatin sample and that all the anomalous behavior occurs right after the quench. The position of the particle trapped in the center of the melted droplet is continuously measured during a succession of several melting and aging. For each measurement the quenching is repeated a few hundred times in order to perform proper ensemble averages.

An example of trajectory obtained with the second set-up is presented in figure 7. When the intensity of the laser is high, the gelatin droplet is in the "sol" phase and the particle fluctuates in an optical trap with a high stiffness. When the intensity of the laser is low, the gelation is occurring and the stiffness of the trap is low (which is the reason why the position fluctuations are bigger). On the red part of the trajectory, the intensity of the laser is high and the gelatin droplet is liquid. On the blue part, the intensity of the laser is low and the gelation is occuring.

Results

In this section we present some results showing that there are no anomalous fluctuations occurring in the ageing of gelatin solution right after a fast quench. We also discuss why this effect that was previously observed is likely to actually be an artefact due to data analysis.

Time evolution of bulk properties and hysteresis

We first performed preliminary measurement of gelatin gelation in bulk (i.e. without the local heating method). We prepare a cell with gelatin solution as described in 2.1. We set the temperature controls of both the objective and of the cell at 37 • C for 30 min to melt the gelatin. We trap one particle and we switch both temperature controls to a given temperature T fb . We wait a few tens of minutes (typically ∼ 30 min) and we set the distance between the particle and the bottom of the cell when there is no more drift due to thermal expansion of the microscope objective. After that, we measure the position of the trapped bead for a long time (e.g. 8 h) at 400 Hz to see the bulk gelation of the sample. For all these measurements, the distance between the bead and the surface is h = 15 µm. The liquid gelatin solution at 37 • C is a Newtonian fluid, and the stiffness k of the trap can be computed directly from the variance of the position x of the bead:

σ 2 x = x 2 = k B T k (4)
For all these measurement the stiffness of the trap was k = 0.46 ± 0.01 pN/µm.

Since the T gel is expected to be around 29 • C, we varied T fb from 31 • C to 27.5 • C. It was found that above 28.5 • C, the gelation does not occur on the time of the experiment and the solution stays liquid, even if its viscosity increases continuously. Below 28.3 • C, the gelation occurs before the end of the experiment. It was estimated that the bulk gelation of the cell volume takes ∼ 260 min at 28.3 • C and ∼ 120 min at 27.5 • C.

Estimating the state ("sol" or "gel") of the gelatin solution, is not trivial, since the fluid can be really viscous without being completely gelled. Qualitatively, the trajectory of the trapped bead starts to be heckled, and the bead sometimes escapes the trapping (see figure 8). An a posteriori test consists in switching off the laser (resulting in switching off the trapping) and letting the sample at T fb for a few more hours (typically over night) to see if the particle slowly fall to the bottom of the cell. If the particle does not fall, the gelatin solution is considered to be fully gelled. To estimate the evolution of the viscosity during the gelation process, we used passive micro-rheology techniques [START_REF] Mizuno | Active and passive microrheology in equilibrium and nonequilibrium systems[END_REF]. The trajectories were divided in portions of ∼ 1 h, and the Power Spectral Density (PSD) was computed for each portion. A long trajectory is required because we need low frequencies to correctly estimate the PSD. We explicitly assume that the bulk aging is slow enough for not perturbing too much the estimation of the PSD when a long trajectory is taken. Or at least, that taking a long time-window will only smooth the rheology result.

As seen in figure 9a, shortly after the decrease of temperature, the PSD is still Lorentzian, as awaited for a particle trapped in a Newtonian fluid at equilibrium [START_REF] Berg-Sørensen | Power spectrum analysis for optical tweezers[END_REF]. The viscous drag coefficient γ = 6πRη (with η the dynamical viscosity of the solution) can be estimated from the value of the cut-off frequency f c = k/(2πγ). Here we find: η = 21 ± 1 × 10 -3 Pa s. As the gelation occurs, the PSD is less and less Lorentzian (see figure 9b), which is the sign that the gelatin solution starts to behave as a viscoelastic fluid [START_REF] Pesce | Microrheology of complex fluids using optical tweezers: a comparison with macrorheological measurements[END_REF].

We plot in figure 10 the evolution of the fitted cut-off frequency at different time after the gelatin solution was set at T fb = 28.5 • C. Even if the spectrum is no longer Lorentzian near the end of the measurement, it seems that the cut-off frequency decreases exponentially. Therefore, the apparent viscosity increase is exponential. From these preliminary measurements, we estimate that the T gel is about 28.3 • C for our gelatin solution at 5 wt%. We chose to work with T fb < 28.3 • C for all the following quenching experiments. As mentioned earlier, gelatin solutions have big hysteretic behaviour [START_REF] Djabourov | Gelation of aqueous gelatin solutions. i. structural investigation[END_REF][START_REF] Djabourov | Gelation of aqueous gelatin solutions. ii. rheology of the sol-gel transition[END_REF]. It follows that the viscoelastic properties of the solution in an important temperature range around T gel cannot be known independently of the sample's history.

Another consequence of the hysteretic behaviour is that the first bulk gelation of the sample must be done in a controlled and reproducible manner. If the sample is let gel for a too long time (generally more than one day), or at a too low temperature ( 22 • C), the first melting/regelling cycles used for the quenching experiment will be different from the following ones (where a "reproducible" state is reached). Especially, in this case the first melting is more difficult to reach. Examples of trajectories are shown in figure 11. One can clearly see some position drifts occurring when the temperature is increased, before reaching a "sol" state where the particle fluctuates in the optical trap. Note that, contrary to figure 7, we use here the first set-up where the trapping laser is different from the heating laser (and the focal point of the two laser are not aligned). Hence, there is nearly no change of the trap stiffness when the intensity of the heating laser is changed. The change of fluctuations amplitude when the heating laser is switched OFF is mostly due to the rapidly increase of gelatin viscosity when gelation occurs.

Difference between ensemble variance and temporal variance in the presence of a drift

We now consider quenching experiment as described in 2.3. We obtain several temporal trajectories of the particles positions for a given quenching temperature T fb < T gel . The important point is to estimate correctly the statistical properties from this set of data. In particular, we are interested in the variance of the position, which has been seen to have an anomalous increase right after the quench [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF].

The correct ensemble variance should be estimated instantaneously at a given time t, by considering the N different trajectories at this time t. If one wants to increase the statistics by taking a small time-window δt, there are at least 3 ways to compute the variance from the set of trajectories. These different ways are schematically represented on figure 12. We call x i (t) the position of the particle for the i th quench at the time t:

• The temporal variance σ 2 time is obtained by estimating the variance over the time δt for each quench, and then averaging over the N quenches:

σ 2 time (t) = 1 N N i=1 1 δt t+δt t (x i (t ) -xi (t)) 2 dt ( 5
)
where xi (t) = 1 δt t+δt t

x i (t ) dt is the temporal mean of x for the i th quench, between t and t + δt.

• The ensemble variance σ 2 ensemble is obtained by estimating the variance over the N quenches at a time t and then averaging over the time-window δt:

σ 2 ensemble (t) = 1 δt t+δt t 1 N -1 N i=1 (x i (t ) -x(t ) ) 2 dt (6)
where

x(t ) = 1 N N i=1 x i (t )
is the ensemble mean of the N trajectories x i (t ) at time t .

• The boxed variance σ 2 box is obtained by taking the N segments of trajectory from x i (t) to x i (t + δt), and then estimating the variance of the whole set of data:

σ 2 box (t) = 1 N δt N i=1 t+δt t (x i (t ) -x (t)) 2 dt (7)
where

x (t) = 1 N δt N i=1 t+δt t
x i (t ) dt is the mean computed on the set of data made of the N segments from x i (t) to x i (t+δt). It is the variance used in references [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF].

Nota Bene: Here to clearly distinguish the role of the time and the ensemble averages we have considered the time as a continuous variable and the number of trajectories as discrete. But experimentally the time is of course also a discrete variable, since we take measurements with a finite sampling frequency. time is computed by estimating the variance of the points in the fuschia box, and then averaging over the trajectories. The ensemble variance σ 2 ensemble is computed by estimating the variance of the points in the green box, and then averaging over the time-window δt. The boxed variance σ 2 box is computed directly by estimating the variance of all the points in the orange box.

Time N trajectories δt

If the system is at equilibrium and δt is big enough to correctly take account of the low-frequency of the motion, all these values should be equal to the equipartition value k B T /k, with k B the Boltzmann constant, T the temperature and k the trap's stiffness.

Unfortunately, when the system is non-stationary (which is the case for an ageing system), these 3 definitions of the variance are not equivalent. Especially, if there's a slow drift existing on each trajectory, the estimations that average over time (i.e. temporal and boxed variances) are likely to show a strong artefact.

To illustrate this effect, we have taken a set of 178 quenches done with the first setup described in 2.2 at 28 • C, sampled at 400 Hz. The parameters were: melting time τ melt = 250 s, melting intensity I melt = 235 W, resting time τ rest = 305 s and trap stiffness3 k = 3.7 pN/µm. One can clearly see on the trajectories that there is a small drift of ∼ 40 nm which occurs right after the quench (see figure 13). Such a drift is often seen for this kind of measurement. We interpret it as a slow relaxation of the gel network, which occurs on a time much smaller than the gelation, but much greater than the heat dissipation. In other words, when the gelation occurs, the particle is trapped in the gel network at a given position. And even if we melt a small droplet, the gelatin network will somehow "remind" this position and pull the particle back to its place when it regels. Here the drift is very visible because the position of the trapping laser is not the same as the position of the locally heating laser. Thus the position where the particle was during the first bulk gelation is not the position where the particle is attracted to when the gelatin is melted. But even when there is only one laser used for both trapping and heating, this drift can occur. It is indeed impossible to verify that the position where the particle is when the sample gelled is exactly the position of the laser, and a drift of only a few nm can be visible. This kind of drift can be avoided by having a more powerful heating laser to completely melt the gelatin on a larger area, as in [START_REF] Gomez-Solano | Nucleation and growth of thermoreversible polymer gels[END_REF]. We then have three characteristic times :

• τ gel the time needed for the gelatin solution to regel completely. It goes from a few hundreds to more than 1000 s depending on the quench temperature T fb .

• τ dynamics the typical time of the particle motion, which is directly 1 fc and evolves from ∼ 5 to ∼ 100 s during the gelation process.

• τ drift the time where the drift is visible, which is typically 1 s for our experiment.

If we take a δt sufficiently small compared to τ drift , the boxed and ensemble variances will give more or less the same result. Whereas, since τ drift < τ dynamics , it is clear that the temporal variance will dramatically underestimate the variance due to the lack of low frequencies signal. Indeed, the temporal variance would require a δt of the order of magnitude of τ dynamics for a correct estimation, which cannot be used because of the drift and the aging. Data are shown on figure 14 for δt = 0.1 s. Now, if the chosen δt is too big compared to the characteristic time of the drift, the boxed variances will start to show an anomalous increase. Data are shown on figure 15 for δt = 1 s. This increase is not a real non-equilibrium effect due to the sol-gel transition, but only an artefact due to data analysis in presence of a slow drift. However, this slow drift is due to the fact that the sample is a gelatin solution, where an elastic network is created in the "gel" phase.

It is nevertheless interesting to see that the correct ensemble variance seems to satisfy the equilibrium equipartition relation at any time after the fast quench, even though there is a clear evolution of the visco-elastic properties with time, and even in the presence of a slow drift at the beginning of the quench:

∀t : σ 2 ensemble (t) = k B T k . ( 8
)
Similar results were seen for different quenches temperatures from 28

• C to 26 • C.
In [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF] it is stated that for longer times after the quench, the variance should decrease because of the elasticity of the gelatin network (as seen figure 1a). This result is not clear. Of course, the particle dynamics becomes arrested after gelation, and at some point the amplitude of its positional fluctuations must decrease in time. At the gel point, the storage modulus and the yield stress of the gelatin sample becomes non-zero, and even for a freely suspended particle the dynamics becomes subdiffusive [START_REF] Larsen | Microrheology of the liquid-solid transition during gelation[END_REF]. Therefore it is clear that the temporal variance calculated for each trajectory should go to zero at long time. However, the ensemble variance calculated instantaneously over several trajectories should not go to zero. Otherwise it would mean that the position where the particle stops to move is the same for every trajectory. Thus, the fact that the variance decreases for long time in previous works is a sign that the ensemble analysis is mixed with some temporal analysis.

It was verified experimentally that for τ rest up to 900 s the correct ensemble variance remains constant. Afterwards, it is not clear that equation 8 is verified, but we do not observe a shrinkage of the distribution of positions.

Correct Estimation of Position Distribution Function

We now want to study not only the variance of the position fluctuations but also their complete Probability Distribution Function (PDF). In order to minimize the risk of slow drifts and to increase the sample frequency, we did new measurements with the second experimental set-up described in 2.2. With this set-up, the trapping and heating laser are the same, and the sample frequency can go up to 10 kHz. However, the calibration of the measured deflection from V to µm requires a supplementary assumption (for example, that the Fluctuation-Dissipation Theorem is verified, as it is done in [START_REF] Fischer | Active-passive calibration of optical tweezers in viscoelastic media[END_REF]). We will start by doing no assumption and plot the results only in arbitrary units.

One must pay attention at which "position fluctuation" is considered, as people often look at the distribution of δx = x -x . The x is the mean of the position x, which can be defined in several ways when the system is not a classic stationary ergodic system. Especially, when one considers a small time-window δt, the correct mean should be the ensemble average x(t) estimated for each time t (as defined in equation 6). But if one takes instead the temporal average xi (t), estimated for each trajectory between t and δt (as defined in equation 5), the results will differ.

As an example, we take the data of 132 quenches at 27.5 • C, sampled at 8 kHz. The parameters are: melting time τ melt = 200 s, melting intensity I melt = 245 W, resting time τ rest = 570 s, and trapping intensity I trap=26 W which corresponds to trap stiffness4 k ∼ 5 pN/µm. We compare the Probability Distribution Function of the positions with a δt = 0.5 s at different times after the quench, when we subtract either the ensemble average (figure 16a) or the temporal average (figure 16b). In the first case, the PDFs are nearly always Gaussian and do not evolve in time. In the second case, the PDFs are always nice gaussians, but with a variance that decreases in time. This is consistent with the previous results showing that the ensemble variance is constant at any time after the quench, whereas the temporal variance decreases logarithmically with the time after the quench. And the variances estimated by doing a Gaussian fit on the PDFs clearly shows the same behaviour (see figure 17). This effect is simple to understand: the trajectories evolve on a time τ gel . This time is much bigger than τ fluc , the typical time of the fluctuations, and δt. On the time-window δt, each portion of trajectory x i (t) can be written x i (t) = xi + δx i (t), where xi is the time average of the i th trajectory over the time-window. When one considers the N trajectory fragments between t and t + δt, the difference between them is mostly due to the averaged value xi of each trajectory fragment, and not to the fast fluctuations δx i (t). Which means that the distribution of all the x i (t) between t and t+δt is nearly the same as the ensemble distribution of the xi . Whereas, the distribution of all the δx i (t) is nothing more than the distribution of the fast temporal fluctuations of one single trajectory.

This difference is very important, as any kind of high-pass filtering (for example a "detrend" function which is often used to suppress slow drifts) done to the trajectories will result in subtracting the temporal average, and thus distort the PDFs estimation.

The experimental results show that the correct estimated PDFs do not evolve in time after a fast quench. Since we have already shown that the correct ensemble variance always verifies the equilibrium equipartition relation, we can conclude that the variance of these PDFs is simply k B T /k. It is again interesting to see that, even if the gelatin solution is aging, its ensemble statistical properties seem to verify relations that are normally verified at equilibrium.

It was also verified with some available data from [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF] that the correct ensemble PDFs are not evolving with time after the quench.

What about heat and Fluctuation Dissipation Theorem?

In previous works [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF] the anomalous fluctuations observed right after the quench were interpreted in terms of heat exchanges between the bath and the particle. Indeed, the heat exchanged between t and t + τ is equal to the variation of the particle's energy ∆U t,τ = ∆U t+τ -∆U t :

Q t,τ = ∆U t,τ = k 2 x 2 (t + τ ) -x 2 (t) . ( 9
)
In particular, the decrease of the variance after the quench has been interpreted as the sign of a heat transfer from the particle to the bath :

Q t,τ = k 2 σ 2 (t + τ ) -σ 2 (t) ≤ 0. ( 10
)
The Probability Distribution Functions (PDF) of the Q t,τ were shown to be asymmetrical for values of t and τ chosen right after the quench (i.e. where the anomalous fluctuations were observed).

A violation of Fluctuation-Dissipation Relation was also observed for times right after the fast quench. It was linked to the non-zero heat exchange by a modification of the Harada-Sasa equality [START_REF] Harada | Equality connecting energy dissipation with a violation of the fluctuation-response relation[END_REF][START_REF] Verley | Fluctuation theorems and inequalities generalizing the second law of thermodynamics out of equilibrium[END_REF] for non-stationary systems:

∞ 1/∆t S x (t, f ) - 2k B T πf Im{ R(t, f )} df = 2| Q t,∆t | k ( 11
)
Where S x (t, f ) is the Power Spectral Density of x and R(t, f ) is the Fourier transform of the linear response function of the position x to a perturbative time-dependent force (these two quantifies are function of the frequency f , but also of the time t since the system is ageing).

All these interpretations comes from the fact that the variance was seen anomalously high right after the quench, and then reduces to the equipartition value after a given time. In particular, the asymmetry and the shape of the PDFs of Q t,τ are simply mathematical consequences of the fact that x(t+τ ) and x(t) have Gaussian PDFs with different variances σ 2 (t + τ ) > σ 2 (t). For example, the exchange Fluctuation Theorem (xFT) that was retrieved with the asymmetrical PDFs of Q is mathematically verified for any random variable defined by y = x 1 -x 2 where x 1 and x 2 are random variables with centered Gaussian distribution of different variances: σ 2

x 1 = σ 2 x 2 . Since we have already shown that, if estimated correctly, the PDFs of x show no anomalous behavior and have a constant variance equal to k B T /k, if follows directly that the PDFs of Q t,τ are symmetrical. Consequently, in average there is no heat exchange between the particle and the bath, for any t and t + τ , and the xFT reduces to a trivial equality because the asymmetry function of a symmetrical distribution is always zero.

Considering the Fluctuation-Dissipation Theorem, one must remind that it is a priori not a good idea to test it in Fourier space. Indeed it is necessary to assume that the system is stationary and ergodic to link the correlation function to the power spectrum with the Wiener-Khinchin theorem [START_REF] Wiener | Generalized harmonic analysis[END_REF][START_REF] Khintchine | Korrelationstheorie der stationären stochastischen prozesse[END_REF]. Therefore, when the system is not stationary, one should in theory look at the proper ensemble correlation function:

EnsCorr xx (t, τ ) = 1 N N i=1 [x i (t) -x(t) ] × [x i (t + τ ) -x(t + τ ) ] (12)
instead of the Power Spectral Density (PSD), which is a temporal quantity. Of course, one can always define a PSD of x i on a given time-window δt for each trajectory S x i (t, f ). And this PSD would be equal to the Fourier Transform (FT) of the temporal correlation of x i computed on the same time-window:

TimeCorr xx (t, τ ) = 1 δt t+δt t [x i (t ) -xi ] × [x i (t + τ ) -xi ] dt = FT{S x i (t, f )}. ( 13
)
But the system needs to be considered stationary and ergodic on the time-window δt, so that the ensemble and temporal correlations should be equal.

Here, the assumption of local stationarity is reasonable since the PSD were computed on 15 s long time-windows (which is short compared to the ∼ 900 s necessary to gel). However, it seems probable that the observed violation of Fluctuation-Dissipation Theorem was only due to the same kind of artifact already responsible for anomalous variance increase (for example: slow drifts for times right after the fast quench), because PSDs are sensible to low-frequency noises. Thus, there is no reason that this apparent violation is linked to an heat exchange, which does not exist anyway. Some ensemble correlation functions of the particle's position are shown in figure 18 for a set of 40 quenches at 27 • C, sampled at 8 kHz. The parameters are: melting time τ melt = 200 s, melting intensity I melt = 270 W, resting time τ rest = 570 s, and trapping intensity I trap=26 W which corresponds to trap stiffness5 k ∼ 5 pN/µm. The data are very noisy, but there is a tendency: the characteristic time increases after the quench, which is reasonable since the gelatin viscosity is also increasing during the gelation. The correlation functions are not simply exponential relaxations, which is consistent with the fact that the PSD are not Lorentzian (as shown in figure 9b).

We also made some experimental tests of Fluctuation Dissipation Theorem (FDT), by looking at the ensemble correlation of x and the response to an Heaviside change of the position of the trap. For these measurements, the position of the trap is changed from X 1 to X 2 at a time t R after the first quench, and the sample is let gel with the particle in X 2 .

Then, for the second quench, the position of the trap is moved back to X 1 at time t R after the quench, and the sample is let gel in X 1 . The procedure is then repeated alternatively. The perturbation introduced by the change of trapping position allows us to compute a normalized response function, averaged over the trajectories:

χ(t R , τ ) = x(t R + τ ) -X initial X final -X initial (14)
where [X initial ; X final ] = [X 1 ; X 2 ] or [X 2 ; X 1 ]. It corresponds to the usual definition of the response function:

χ(t) = x(t) perturbed -x(t) unperturbed perturbation amplitude . ( 15
)
We use X final -X initial which is proportional to the perturbation amplitude. And we simply take X initial as the average value of the unperturbed trajectory, because the mean position of the bead is constant and equal to the position of the trap if there is no perturbation 6 . If the FDT is verified, the response function should verify:

χ(t R , τ ) = 1 - k k B T EnsCorr xx (t R , τ ) (16)
Some data are presented in figure 19 for 50 quenches at 26 • C, sampled at 8 kHz. The parameters are: melting time τ melt = 200 s, melting intensity I melt = 270 W, resting time τ rest = 570 s, and trapping intensity I trap=26 W which corresponds to trap stiffness7 k ∼ 5 pN/µm. The values of X 1 and X 2 are estimated by computing the mean position of the bead when the gelatin is melted (which gives alternatively X 1 and X 2 ). The exact value of k B T /k was extracted from the variance of the position PDFs computed before changing the position of the trap. These measurements are a bit noisy because it requires a lot of statistics to compute a proper ensemble correlation function, but no apparent violation of the FDT was found for the times tested. We didn't test the Fluctuation Dissipation Theorem for times t R taken shortly after the quench, because the ensemble correlation shows a characteristic time which is very short at this time (see figure 18). It is then more difficult to compute a proper ensemble correlation right after the quench, than when the viscosity of gelatin has already started to increase. We also didn't compute the response function by varying t R for a fixed t R + τ , because it would require a lot of time to do the experiments. Indeed, each set of t R requires one day of measurement to compute χ(t R , τ ), and the sample cannot be kept a lot of days without degrading.

Conclusion

In conclusion, we have locally studied the gel transition of a gelatin solution. We were unable to reproduce the results of previous works [START_REF] Gomez-Solano | Nonequilibrium fluctuations of a Brownian particle[END_REF][START_REF] Gomez-Solano | Heat fluctuations in a nonequilibrium bath[END_REF][START_REF] Gomez-Solano | Fluctuations, linear response and heat flux of an aging system[END_REF], but we have identified some experimental and data analysis artifacts which may explain the effects previously observed. In particular we have analyzed the effect of time-windows on proper ensemble averages, which are important to study aging systems.

We have shown that in the hysteresis range of temperature (28.3 • C < T < 36 • C), bulk gelation can occur on very long times, and viscoelastic properties gradually appear. The characteristic time of the particle trapped in the bulk-gelling sample was seen to decrease exponentially before the gelation (whereas the viscosity evolves logarithmically after the gelation).

For fast quenches of a small droplet of gelatin solution, we have found that the Probability Distribution Functions of the position of the trapped particle do not evolve with time after the quench, even if the gelatin sample is undergoing aging and the viscoelastic properties are clearly evolving. Moreover, these PDFs show equilibrium-like properties, being Gaussian with a variance equal to the equipartition value k B T /k. These results seem not so surprising a posteriori, since it was already observed in the previous works that, after ∼ 15 s the Brownian motion of the trapped particle behaves like in equilibrium with the thermal motion of the gelatin chains. Only the very first seconds after the quench showed anomalous behavior, which was strange, because the complete gelation occurs on much larger scales (∼ 900 s). It however remains striking that a system which is clearly not stationary because of aging has ensemble properties which are stationary.

For systems which are not ergodic or not stationary, time properties can be very different from ensemble properties. And it was also shown that some artifacts (like slow drifts) or analysis bias (like high-pass filter) can greatly modify the results if ensemble properties are estimated on time-windows. Therefore, one must be very careful when studying statistical properties of an aging system. This kind of problems had already arisen for other aging systems. For example, it was already shown in [START_REF] Jop | Experimental study of out-of-equilibrium fluctuations in a colloidal suspension of laponite using optical traps[END_REF] that increase in effective temperature previously seen in suspension of Laponite [START_REF] Greinert | Measurement of effective temperatures in an aging colloidal glass[END_REF] were in fact artifacts due to analysis methods.

Finally, in agreement with the absence of anomalous position fluctuation after the fast quench, no heat exchange, nor violation of the Fluctuation Dissipation Theorem was seen, as it would be expected in an equilibrium medium.

Figure 1 :

 1 Figure 1: (a) Evolution of the normalised variance kσ 2 x /k B T of the position fluctuations of one particle trapped in gelatin solution (10 wt%) or glycreol, quenched at 26 • C, for different times t after the quench. (b) Evolution of the Probability Distribution Function of the position fluctuations of the particle trapped in gelatin solution.

Figure 2 :

 2 Figure2: Probability Density Function of the normalised heat q exchanged during τ = 30 s computed at different times t after the quench in gelatin.

Figure 3 :

 3 Figure 3: Passive Power Spectral Densities of the position fluctuations (color points) and Fourier transform of the active response function (black dashed-lines) computed at different times after the quench. (a) For 0 s < t < 15 s. (b) For 30 s < t < 45 s. (c) For 75 s < t < 90 s. (d) For 1200 s < t < 1215 s. If the Fluctuations Dissipation Theorem is verified, the two quantities should be equal, which is not the case for low-frequency in (a).

Figure 4 :

 4 Figure 4: Schematic representation of the cell used to trap particles in gelatin solution (view from the side). The microscope slide is ITO-coated, which enables us to heat the cell by sending an electrical current through the glass surface.

Figure 5 :

 5 Figure5: Schematic representation of optical tweezers set-up used to trap two particles in gelatin solution with a green laser. An infrared laser diode is used to melt the gelatin. "M" is a mirror and "DM" are dichroic mirrors.

Figure 6 :

 6 Figure6: Schematic representation of optical tweezers set-up to trap one particle in gelatin solution. An infrared laser diode is used to trap and to melt the gelatin. The deflection of a He-Ne laser induced by the trapped bead is measured by the position sensing diode (PSD). The white light source and camera are only used for direct visualization but not for the measurements. "BS" is a beam-splitter, "IF" and interferential filter to suppress the infrared beam and "DM" are dichroic mirrors.

Figure 7 :

 7 Figure7: Trajectory of one particle trapped in a gelatin sample kept at T fb = 27.5 • C, around the quench. On the red part of the trajectory, the intensity of the laser is high and the gelatin droplet is liquid. On the blue part, the intensity of the laser is low and the gelation is occuring.

Figure 8 :

 8 Figure 8: Evolution of the position of one trapped particle, in gelatin solution kept at 28.3 • C (after being melt at 37 • C for 20 min). At the end of the trajectory, gelation occurs and the particle is moved away from the optical trap.

  70 min after the temperature change.

  320 min after the temperature change.

Figure 9 :Figure 10 :

 910 Figure9: Power Spectral Density of the position of one particle trapped in gelatin solution kept at 28.5 • C (after being melt at 37 • C for 30 min). The PSD is estimated over a timewindow of 1 h. Shortly after switching the temperature the gelatin solution is still a Newtonian fluid and the PSD is Lorentzian. After some time, visco-elastic effects appear and the PSD is no longer Lorentzian.

  3 h after starting the cycles.

Figure 11 :

 11 Figure11: Examples of trajectories for melting/regelling cycles with a gelatine sample gelled at T room = 24 • C for a long time. The red dashed-lines indicate when the heating laser is switched ON, and the blue ones when it is switched OFF. At the beginning, the melting is more difficult to reach, after a given time, the cycles look "reproducible".

Figure 12 :

 12 Figure 12: Schematic representation of the different ways to estimate the variance for a set of N trajectories with a time-window δt. The temporal variance σ 2time is computed by estimating the variance of the points in the fuschia box, and then averaging over the trajectories. The ensemble variance σ 2 ensemble is computed by estimating the variance of the points in the green box, and then averaging over the time-window δt. The boxed variance σ2 box is computed directly by estimating the variance of all the points in the orange box.

Figure 13 :

 13 Figure 13: 20 first trajectories for a quench at T fb = 28 • C, sampled at 400 Hz. A slow drift of ∼ 40 nm is clearly visible during the first ∼ 1 s. After that, the position only oscillates randomly around a mean value.

  Ensemble and boxed variances.

Figure 14 :

 14 Figure 14: Different variances computed for δt = 0.1 s and normalised by the equilibrium value k B T /k. The ensemble and boxed values are nearly equal and seem to be close to the equilibrium value at any time after the quench. Whereas the temporal value is clearly below the equilibrium value and decreases logarithmically with time after the quench.

2 Figure 15 :

 215 Figure15: Ensemble and boxed variances computed for a δt = 1 s and normalised by the equilibrium value k B T /k. The boxed variance clearly shows an anomalous increase at small times after the quench, which is the effect reported in previous works.

  When one subtracts the temporal average.

Figure 16 :

 16 Figure 16: Evolution of the Probability Distribution Function of the position fluctuation δx depending on the definition taken for the subtracted average. The PDFs are computed on a time-window δt = 0.5 s for different times after the quench, going from t = 0 s (blue curves) to t = 540 s (red curves).

Figure 17 :

 17 Figure17: Evolution of the variance estimated by fitting the PDFs with a Gaussian, at different times after the fast quench. When subtracting the correct ensemble average the variance is constant. When subtracting the temporal average, the variance decreases almost logarithmically with the time after the quench.

Figure 18 :

 18 Figure 18: Normalised ensemble correlation function for a quench at 27 • C. Here we keep t fixed and we vary τ from -10 s to 0 s. The normalisation is done by dividing EnsCorr xx (t, τ ) by the value of k B T /k extracted from the variance of the position PDFs.

Figure 19 :

 19 Figure 19: Normalised response function χ(t R , τ ) and ensemble correlation function for t R = 100 s after the quench, and τ going from 0 to 10 s.

This is the power measured on the beam before the microscope objective, so the "real" power at the focal point should be smaller, due to the loss in the objective.

Actually, in the case where the same laser diode is used for trapping and melting the droplet, the temperature of the quench is a little bit above T fb because of the absorption of the laser. Since the power of the laser is low, this increase is less than 1 • C and can easily be compensated by lowering T fb accordingly.

The trap stiffness is measured when the gelatin sample is completely melt and kept at constant temperature T = 37 • C, before the first bulk gelation.

The trap stiffness was measured in water (where viscosity is known) for the same laser intensity.

The trap stiffness was measured in water (where viscosity is known) for the same laser intensity.

One could also take x(t R ) to guarantee that χ(t R , 0) = 0, but it wasn't necessary here.

The trap stiffness was measured in water (where viscosity is known) for the same laser intensity.