Lauréline Seurre
email: laureline.seurre@etu.uphf.fr

Sofiane Ghenna

Héléne Arena

Caroline Soyer

Sébastien Grondel

Cedric Plesse

Giao Tran

Minh Nguyen

Frederic Vidal

Eric Cattan

H Aréna

G T M Nguyen

Memory effect and dynamics in PEDOT:PSS-based actuators under DC voltage

Keywords: Electroactive polymer, electronically conducting polymer, memory effect, relaxation, DC voltage

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

There is a growing interest in future soft micro-robotics, which requires the use of highly compliant, flexible, and good adaptable materials A category of these micro-robots will be used as micro-tools to improve the precision of the surgeons' work or to explore the human body, provided the tools are contained within volumes of a few cubic millimeters. The development of new generations of flexible materials that can be used as micro actuators and micro sensors is therefore necessary. Ionic electroactive polymers are especially attractive owing to their striking characteristics, vital for medical applications, of low driving voltages and operation under ambient conditions. Various types of ionic electroactive polymers operating under low voltage conditions have been developed, including conducting polymers [START_REF] Baughman | Conducting polymer artificial muscles[END_REF][START_REF] Kaneto | Artificial muscle': Electromechanical actuators using polyaniline films[END_REF][START_REF] Smela | Conjugated polymer actuators for biomedical applications[END_REF][START_REF] Plesse | Conducting polymer artificial muscle fibres: Toward an open air linear actuation[END_REF][START_REF] Vidal | Long-life air working conducting semi-IPN/ionic liquid based actuator[END_REF][START_REF] Khaldi | Conducting interpenetrating polymer network sized to fabricate microactuators[END_REF][START_REF] Khaldi | Smarter Actuator Design with Complementary and Synergetic Functions[END_REF][START_REF] Maziz | Top-down Approach for the Direct Synthesis, Patterning, and Operation of Artificial Micromuscles on Flexible Substrates[END_REF] , ionic polymer-metal composites [START_REF] Asaka | Bending of Polyelectrolyte Membrane-Platinum Composites by Electric Stimuli I. Response Characteristics to Various Waveforms[END_REF][START_REF] Nemat-Nasser | Electromechanical response of ionic polymer-metal composites[END_REF][START_REF] Tiwari | The state of understanding of ionic polymer metal composite architecture: A review[END_REF] , bulky gel polymers [START_REF] Fukushima | Fully plastic actuator through layer-by-layer casting with ionicliquid-based bucky gel[END_REF][START_REF] Lee | Bucky gels' for tailoring electroactive materials and devices: The composites of carbon materials with ionic liquids[END_REF][START_REF] Asaka | Highly conductive Sheets from Millimeter-Long Single-Walled carbon nanotubes and ionic liquids: Application to Fast-Moving, Low-Voltage electromechanical actuators cperable in airs[END_REF][START_REF] Imaizumi | Driving mechanisms of ionic polymer actuators having electric double layer capacitor structures[END_REF] . However, very few of them were developed to be subsequently miniaturized using microelectronic fabrication technologies such as photolithography patterning, thin films deposition and plasma etching, with thicknesses less than a few microns equivalent to those of the rigid micro-actuators traditionally used for microsystems production [START_REF] Smela | Controlled folding of micrometer-size structures[END_REF][START_REF] Jager | Microfabricating conjugated polymer actuators[END_REF] . Soft microsystems, including high-performance thin flexible actuators and sensors, and integrated proximity electronics on a flexible base, are essential elements for future soft micro-robots. Very thin conducting polymer poly (3,4-ethylenedioxythiopene) (PEDOT) and more recently poly(3,4-ethylene dioxythiophene):poly(styrene sulfonate) (PEDOT:PSS) based-actuators have been integrated into soft micro-structures with gold remote contacts [START_REF] Rohtlaid | PSS-based micromuscles and microsensors fully integrated in flexible chips[END_REF][START_REF] Seurre | Demonstrating Full Integration Process for Electroactive Polymer Microtransducers to Realize Soft Microchips[END_REF] . They have also received significant attention because they exhibit large strains and forces with a driving voltage of less than 2 V.

To compare the different types of actuator devices, their electromechanical performances, among other things, have primarily been studied with a square alternating current. However, this is not usually the operating mode of the actuators involved in micro-robot applications. A direct current (DC) signal or a voltage set-point step is often tuned to obtain the desired position of the actuator. Their motions are not periodic and do not alternate over time, but they should remain in the same position for several seconds or minutes. Applications, for which a quasistatic voltage is necessary, include the displacement of a micro-camera on the tip of an endoscope, the actuation of a cochlear implant during surgery, or closing micro-tweezers for handling human tissue. These applications show the great potential of conducting polymer basedactuators, and understanding their working mechanisms under DC voltage is becoming increasingly important for developing and controlling microdevices. In this way, knowing the dynamics of the actuator under DC voltage to reach a strain difference (SD) or force setpoint is necessary to understand the phenomenon. This solicitation results in the onset of a residual strain difference (RSD) when the supply voltage is shut down, i.e. when the capacitor (actuator) is fully discharged. The dynamics under DC voltage depend heavily on the materials used in single, bi, or trilayers and on whether the measurements are made in solution or air. Extensive studies have been carried out on a polypyrrole matrix, in particular by Otero et al. on the current, charge, and strain dynamics after previously subjecting the polymer film to cathodic potentials for long periods of time. The interpretation given is that the polymeric structure can be compacted either by applying high cathodic potentials or by maintaining those potentials for long periods of time. In this sense, a good combination of potential and time allows the degree of structure closure to be improved without any decrease in electrical properties due to degradation processes [START_REF] Grande | Conformational movements explain logarithmic relaxation in conducting polymers[END_REF][START_REF] Grande | Intrinsic asymmetry, hysteresis, and conformational relaxation during redox switching in polypyrrole: A coulovoltametric study[END_REF][START_REF] Otero | Reinterpretation of polypyrrole electrochemistry after consideration of conformational relaxation processes[END_REF][START_REF] Sendai | Anisotropic strain and memory effect in electrochemomechanical strain of polypyrrole films under high tensile stresses[END_REF] . Sendai et al. [START_REF] Sendai | Anisotropic strain and memory effect in electrochemomechanical strain of polypyrrole films under high tensile stresses[END_REF] studied electrochemomechanical strain under tensile stress in polypyrrole films, and discussed the results in terms of the relaxation of the anisotropic strain and the memory effect introduced during creeping. The actuator dynamics is affected by memory effect described also in [START_REF] Yamaura | Memory effect of electrical conductivity upon the counteranion exchange of polypyrrole films[END_REF][START_REF] Vorotyntsev | Memory effects in functionalized conducting polymer films: Titanocene derivatized polypyrrole in contact with THF solutions[END_REF] for polypyrrole matrix with pendant titanocene dichloride centers, by applying a combination of different experimental conditions. It has been established that both electroactive components of the film, the polymer matrix and the immobilized centers, contribute to the memory effects. The studies on dynamics with ionic polymer metal composites (IPMCs) trilayer actuators [START_REF] Porfiri | An alternative explanation of back-relaxation in ionic polymer metal composites[END_REF][START_REF] Vunder | Mechanical interpretation of back-relaxation of ionic electroactive polymer actuators[END_REF][START_REF] Liu | Ion transport and storage of ionic liquids in ionic polymer conductor network composites[END_REF] focus more on the back-relaxation effect: actuator under DC voltage that relaxes slowly back towards its initial shape instead of holding its bent state. In this case, maintaining a static position other than neutral due to hysteresis and creep is difficult [START_REF] Jo | Recent advances in ionic polymer-metal composite actuators and their modeling and applications[END_REF] , hindering their precise control.

To our knowledge, there is very little applicable data to map the behavioral dynamics of conducting polymer-based actuators when a DC voltage is applied for several minutes. The objective of this paper is to identify the time-dependent correspondence between the transient input DC voltage and the shape of the conducting polymer-based actuators, and propose some interpretations. First, the innovative method used to fabricate the extremely small trilayer actuators tested (with a thickness of 35 µm and a surface area of a few square millimeters) is explained, and their traditionally accepted mode of operation (redox process) is quickly reviewed. Then, the specific measuring protocol related to the need to apply electrical voltages successively with the same polarity is briefly described. No back-relaxation effect on the SD is observed for periods as long as 10 min under voltage. However, in this paper, the study focuses more on the dynamics of relaxation under a short-circuit. The RSD is described, the evolution of SD and RSD following a repeated applied voltage is analyzed, and a macroscopic explanation of the phenomenon is proposed.

MATERIAL AND METHODS

Material

Poly(ethylene glycol) methyl ether methacrylate (PEGM, Mn = 500 g.mol-1), poly(ethylene glycol) dimethacrylate (PEGDM, Mn = 750 g.mol -1 ), ammonium persulfate (APS, 98 %) and cyclohexanone (>99.8 %) were obtained from Sigma Aldrich and used as received. Poly(3,4-ethylenedioxythiophene) polystyrene sulfonate (PEDOT:PSS) aqueous solution (Clevios PH1000, solid content 1.0 -1.3 wt%) was purchased from Heraeus Precious Metals GmbH & Co. Nitrilebutadiene rubber (NBR), initiator dicyclohexyl peroxydicarbonate (DCPD), and 1-ethyl-3-methylimidazolium bis(trifluoromethanesulfonyl)imide (EMImTFSI 99.9%) were used as supplied from LANXESS, Groupe Arnaud and Solvionic, respectively.

Micro-actuators fabrication

The trilayer actuators composed of a PEO-NBR ion storage membrane (ISM) sandwiched between two electroactive PEDOT:PSS-PEO layers were fabricated using the layer stacking method [START_REF] Seurre | Demonstrating Full Integration Process for Electroactive Polymer Microtransducers to Realize Soft Microchips[END_REF] . The PEDOT:PSS-PEO casting solutions were obtained by mixing PEO (poly(ethylene oxide)) precursors (mPEG, methacrylic functionalized PEG) composed of 50 wt% PEGM and 50 wt% PEGDM as the monomer and the crosslinker, respectively, with an aqueous dispersion of PEDOT:PSS (40wt% with respect to the final electrode, i.e. 60wt% PEDOT:PSS, equivalent to 1.15 wt% solid content). APS was also added to the solution (3 wt% with respect to mPEG) as a radical initiator for the precursors. The solution was stirred until dissolution. The resulting solution was cast (0.05 ml.cm -2 ) onto a glass slide and placed on a heating plate at 50°C, to evaporate the water. The ion storage membrane located between the two polymer electrodes is based on a semiinterpenetrating polymer network layer composed of a PEO network (50 wt%) and a linear NBR (50wt%). First, the NBR solution was prepared by dissolving NBR in cyclohexanone to obtain a concentration of 20 wt%. The PEO precursors, consisting of PEGDM (25 wt% with respect to the PEO network) and PEGM (75 wt% with respect to the PEO network) were added to the NBR solution and stirred for 30 min. The radical initiator DCPD, 3 wt% with respect to the PEO network, was then added to the solution. The final solution was stirred until complete homogenization, then degassed. During the next step, the reactive mixture was spin coated (3000 rpm -1000 rpm/s -30 s) onto the first PEDOT:PSS-PEO electrode layer and pre-polymerized in a closed annealing chamber under continuous nitrogen flow for 45 min at 50°C, to initiate the formation of the PEO network. The second PEDOT:PSS-PEO electrode was fabricated on top of the PEDOT:PSS/PEO-NBR/PEO bilayer in the same way as the first electrode: new solutions were prepared, cast and solidified at 50°C by evaporating the water. The resulting trilayer actuators were then placed in a closed annealing chamber and the final heat treatment was carried out at 50°C for 3 h followed by post-curing at 80°C for 1h under a continuous nitrogen flow. The thicknesses of each PEDOT:PSS/PEO electrode and the NBR/PEO layer obtained were 6µm and 13 µm, respectively. The trilayer structures were then micro-patterned as microbeams using laser cutting and measured between 5 mm and 6 mm in length (L), 1 mm and 1.43 mm wide (w), and 35 ± 1 µm thick (h). The micro-actuator beams were then swollen in ionic liquid (IL) EMImTFSI for 72 h to incorporate the ions necessary for the redox process and to perform characterizations at a later stage. To simplify reading afterwards, the "micro-actuators" will be called "actuators".

RESULTS AND DISCUSSION

Micro-actuator operation

Electronically conducting polymer-based actuators are able to store energy; therefore, they can be assimilated to capacitors. For a better understanding of the following discussion, this part reminds the operation work of these actuators. The principle of actuator operation is based on an oxidation-reduction process of the conductive polymer [START_REF] Otero | Biomimetic electrochemistry from conducting polymers. A review[END_REF][START_REF] Otero | Biomimetic conducting polymers: Synthesis, materials, properties, functions, and devices[END_REF][32][START_REF] Pei | Electrochemical applications of the bending beam method; a novel way to study ion transport in electroactive polymers[END_REF][START_REF] Otero | Electrochemomechanical properties from a bilayer: polypyrrole / non-conducting and flexible material -artificial muscle[END_REF] , going from a doped state (removal of electrons, which corresponds to the oxidized state), to an undoped state (insertion of ions which corresponds to the reduced state). The two conducting polymer layers successively undergo doping/undoping when they are electrically powered up. This electrochemical process, in both conducting polymer layers, leads to structural changes 35- 37 . In our case, the reduction of the p-doped conducting polymer lead to the dissociation of the polymer backbone (PEDOT 0 ) with its counterion, PSS -. PSS -is retained in the matrix because of its large size, which hinders its movement [START_REF] Löffler | Organic bioelectronic tools for biomedical applications[END_REF][START_REF] Park | Redox-active charge carriers of conducting polymers as a tuner of conductivity and its potential window[END_REF] . Therefore, charge balancing is mainly achieved by the displacement of cations from the electrolyte (EMITFSI) in the PEDOT:PSS/PEO electrodes and in the NBR/PEO electrolyte support membrane between the two electrodes (Figure 1). The disappearance of positive charges on the PEDOT backbone makes the cations come into the cathode, while the appearance of positive charges on the PEDOT backbone at the anode makes the cations leave. When using IL, the electrochemical reactions leading to the expansion of one electrode during the reduction and simultaneously to the contraction of the other electrode during the oxidation are described by equation (1), and are responsible for the bending of the actuator:

Oxidation (PEDOT + :PSS -) + e -+ EMI + (PEDOT 0 , PSS -:EMI + ) (1)
Oxidized state Reduction Neutral state

The EMI + are the mobile cations of the electrolyte governing actuation [START_REF] David | Applications of ionic liquids in polymer science and technology[END_REF] , and the TFSI -are the immobile anions of the electrolyte, trapped between the chains of the PEDOT skeleton in the oxidized state, and therefore not participating in the reaction

Set-up and methods for characterization

To perform strain measurements, the actuator beams were hold between two glass plates covered with copper wires to apply the DC voltage. The maximum voltage was 1.75 V to avoid trilayer damage and it was decreased in 0.25 V increments from 1.75 V to 0.5 V. The experiments were performed from the side, in air, using a Dino Lite micro camera. The strain was recorded by measuring the curvature radius r of the beam on photos taken regularly. The radius length was obtained using three points along the edge of the beam. The results were obtained by performing an electrical applied voltage for 10 min at each step. This duration was identified as allowing a level strain to be reached in all cases, i.e. the strain did not change by more than 0.1%/min of SD after 10 min to 20 min of applied voltage. After being switched on for 10 min, each actuator was short-circuited to discharge the charge accumulated by the capacitor, resulting in a peak current discharge lasting for a few seconds. The next measurement was carried out after 10 min of rest to give the beam time to return to its initial position. This time was defined as sufficient, as beyond this the relaxation strain reached a speed of 0.005%SD/min. Unless otherwise specified, this protocol was applied for the following measurements.

The strain difference (SD) is defined by h/r, where h is the thickness of the actuator. This study is focused on the actuator SD after applying an electrical voltage but also after switching off the electric voltage (short circuit). A particularity of the measurements presented is that the polarity remained unchanged when switching from one electrical voltage to another, the actuator is never subjected to an electrical voltage of reverse polarity. This situation is intended to simulate, for example, an actuator placed in a cochlear implant that must be deformed in a single direction to bend the implant in the cochlea as much as possible.

Dynamics of actuator strain under DC voltage

Figure 2b shows an untouched PEDOT:PSS-based actuator, before any electrical applied voltage. First actuation was achieved directly at 1.75 V and maintained for 10 min and a picture was taken to determine the actuator SD. During the first actuation, the PEDOT charges disappear at the cathode, making the cations of the IL move toward the cathode (Figure 1), while at the anode the appearance of positive charges on the PEDOT backbone makes the cations leave this electrode. It creates, in this case, an SD of 2.0% (Figure 2c). The conductive polymer electrodes were then immediately shortcircuited. A new picture was taken after 10 min (Figure 2d) and a residual strain difference (RSD) of 0.6% was determined from the new radius of curvature, characterizing a beam that does not return to its initial vertical position. The actuator finds itself in unstable equilibrium position as, despite everything, it continues to move very slowly after the 10 min shortcircuit as it will be seen later (Figure 3b). Immediately after this first measurement, a second measurement was carried out in the same way (Figure 2e&f). However, when this new measurement at 1.75 V began, the actuator was in the state of strain in Figure 2d (0.6%). The obtained total strain was 2.4% but the strain obtained between its initial and final position was 1.8%, which is close to the value found during the first measurement (2.0%). Therefore, the actuator produces similar SDs under a same electrical voltage, but its final position in space is not the same. A new picture was taken after 10 min of short-circuit (Figure 2f) and a second residual strain difference (RSD2) of 0.97% was determined The non-zero RSD (Figure 2d&f) indicates that the actuator memorizes part of the strain obtained under the previous voltage and that this memory effect is dependent on the accumulation of energy supplied to the actuator in the same direction of operation. For an actuator that has not been powered up, it can be assumed that the ions of the IL are homogeneously distributed in the actuator. During the first applied voltage, electrons are removed from the PEDOT:PSS on the anode side (oxidation) and injected into the cathode side (reduction). The cations near the interface and in the ISM gradually move toward the ECP (the cathode), to compensate the disappearing of the positive charges on the reduced PEDOT backbone, while the cations in the anode are expelled toward the ISM due to the appearance of positive charges on the doped PEDOT backbone (Figure 1). The cations are then concentrated in or near the cathode. This concentration generates an increase in its volume, which causes spacing between the polymer chains of the ECP at the cathode, and the leaving of the cations at the anode generates its contraction, and thus the curvature of the whole actuator (Figure 2c&e). The cations are then subjected to a certain pressure that keep them concentrated in or near the cathode.

To better understand the influence of applied voltage on the RSD, the evolution of the RSD was recorded over 25 h, starting with short-circuit performed after a 10 min of actuation at 1.75 V (Figure 3b). This experiment provided us with the strain dynamics over a period longer than 10 min, which was compared to evolution of the RSD over time for a sample that had not been subjected to any voltage, having simply been mechanically deformed (Figure 3c). To obtain this mechanical strain, the end of the actuator was bent around a cylindrical rod, to reproduce a strain similar to that obtained under electrical tension. This strain was maintained for 10 min, after which the locking system was removed and successive images were recorded for 35 min (Figure 3c). In this case, the RSD decreased sharply and reached zero in approximately 2000 s, indicating that the actuator had returned to its initial vertical position. Its relaxation dynamics can be described by three exponential decay functions, which correspond to a three phases relaxation: a first extremely fast phase with a constant time of only 0.3 s, a second slower phase with a constant time of 5 s, and a final even slower phase with a constant time of 524 s. 82% of the SD was lost in 13 s, and the beam returned to its initial vertical position after 30 min.

In the case of a relaxation after an applied voltage of 1.75 V (Figure 3b), even after 25 h the actuator had not returned to its vertical position. An unstable state was established and the actuator only returned to its initial vertical position after several days. Its relaxation dynamics followed a two decreasing exponential functions, corresponding to a fast relaxation with a constant time of 1 = 13 s and a slower relaxation with a constant time of 2 = 7937 s. The time required to lose 50% of the strain was 26 s, and in 4 h of relaxation the RSD reached 0.12%.

When a short-circuit is caused (actuator switching from the state in Figure 2c to the state in Figure 2d), a combination of at least two physical phenomena can be considered: (i) the reflux of the IL cations, and (ii) the purely mechanical relaxation of the beam.

i.

The short circuit enable the electrons in the ECP electrodes to flow out resulting in a current discharge peak. This charge movement causes the electrons to recombine on the polymer chains releasing the cations concentrated in or near the cathode [START_REF] Löffler | Organic bioelectronic tools for biomedical applications[END_REF] . The cations flow back toward the ISM and the anode, changing the state of contraction and expansion of the two ECP electrodes and causing the beam to relax (Figure 2d&f, Figure 3b). As the cations are within a relatively dense matrix of polymer chains, they are located in an interstitial position between the polymer chains. Therefore, they are slightly trapped and their limited mobility slows down the reflux. Consequently, as long as a reverse electrical voltage is not applied, the cations do not have enough energy to flow back quickly.

ii.

When a short circuit is caused, the beam releases the stored mechanical energy (viscoelasticity of the material).

The ECP electrodes extended on one side and compressed on the other side, try to return to their initial position naturally by participating in the expulsion of the IL located between the polymer chains of the expanded cathode. The comparison of the curves in Figure 3b&c shows that the reflux of the IL ions slows the return of the actuator to the initial position: the constant time in the first phase of purely mechanical relaxation is only 0.3 s, i.e. 47 times faster than during relaxation following an electrical applied voltage and a short-circuit. In the last relaxation phase, the return to the initial position is slower in both situations, but purely mechanical relaxation is still 15 times faster. Therefore, in the first phase of relaxation (Figure 3b), the viscoelasticity of the material is the main contributor to the return of the microbeam to its initial position. This effect is then strongly slowed down by the IL reflux and then an unstable balance is established between the two phenomena. Only the presence of a high concentration of IL near and/or in the cathode can explain the maintenance of a curvature opposing the elastic strain of the beam, thus preserving part of the strain caused by the voltage applied to the beam.

When the actuator was subjected to a second applied voltage at 1.75 V, immediately after the previous 10 min short-circuit (Figure 2e), the total strain achieved relative to an initial vertical position was 2.4%. However, the final strain achieved takes into account the initial state of strain of the beam (RSD1 of 0.6% after the first measurement Figure 2d). When it is subtracted from the total strain, the SD becomes 1.8% (Figure 2). Therefore, a loss of 0.2% (2% to 1.8%) was observed after the second applied voltage. After the 10 min short-circuit, the RSD2 obtained was 0.97% (Figure 2f), which corresponds to a 0.37% increase in RSD compared to the previous RSD1. This repeated actuation (2 × 10 min) resulted in the accumulation of ions in the cathode. Even though the actuator still had a margin of strain, it had to fight against the difficulty of deforming the polymer (maximum elongation of the polymer chains) during the second actuation. This could explain why with the same energy, the SD obtained for the second actuation was slightly lower (by 0.2%). This problem of the mechanical strain limit of the structure is accentuated by the fact that during this second strain the actuator starts from an initial non-zero strain (0.6%, Figure 2d). Moreover, during the first application of an electric voltage, the mobility of the cations was better because they were moved over a "large" distance (from a state where the cations were homogeneously distributed in the actuator to a state where the cations were concentrated in the cathode). With the second applied voltage, the efficiency was lower because the cations had more limited mobility as they were already at the ISM/ECP interface or in the ECP. Conversely, the RSD was higher for the second measurement than for the first one indicating that the ion reflux was less efficient and that the concentration of IL remaining in or near the cathode was higher, thus maintaining the polymer chains in a swollen state at the cathode. This second applied voltage allowed more IL to be blocked between the polymer chains on the cathode side. Therefore, the unstable state of equilibrium, mentioned above is displaced spatially.

These results show that the total strain and the residual strain of the beam, under the applied voltage and time conditions described above, produced a positioning of the beam dependent on the prior distribution of the ions inside the actuator. Beyond a few minutes of actuation, the backflow of the IL following the short-circuit, disturbs the pure elastic mechanical return of the actuator, keeping it in a quasi-stable position and partially memorizing the action performed previously.

Strain difference and residual strain difference under various DC voltage

In the experiments described previously, the only electrical voltage applied was 1.75 V for periods of 10 min to deform the actuator. If such materials were used to actuate micro-robots, a succession of electrical solicitations with various voltage amplitudes and a return to zero between each measurement should be possible. This situation has therefore been taken into account to some extent in the measurements presented in Figure 4. A fresh actuator was used, that had never been electrically energized. Therefore, the initial distribution of ions inside this material was assumed to be homogeneous. The experiment presented in Figure 4b started with a voltage of 1.75 V for 10 min, followed by a short circuit for 10 min, and this protocol was repeated by gradually decreasing the electrical voltage by 0.25 V, but always maintaining the same polarity (see chrono-voltamogram). The SD was calculated as described above, taking into account the RSD before actuation. A second series of SD measurements was performed on the same sample under the same conditions but after about 20 h of rest for the sample. This rest allowed time for part of the IL to be redistributed in the material and to reduce the RSD to around 0.1% (Figure 3b). For the two measurements (first and second series) conducted at a voltage of 1.75 V, it can be seen that the SD values differed by only 0.3% (Figure 4b) and are not on the axis as the line formed by the SD obtained for the other voltages (1.5 V to 0.5 V). These two SDs at 1.75 V, obtained for the initial states where the distribution of the IL was particularly homogeneous are different from the other SD values obtained at the lower electrical voltages (1.5 to 0.5 V). This result confirms the hypothesis described above, which indicates that if the cations distribution is initially homogeneous or has regained a certain homogeneity in the actuator due to resting the material, then the strain is amplified. The offset of 0.3% between the two SDs at 1.75V confirms the results in Figure 2b&d (SD1b-SD1d = 0.2%).

In our opinion, these effects are closely related to the dynamics of cation movement, which is favored when they are homogeneously or relatively homogeneously distributed in the material. Conversely, this dynamic is more limited when the cations are already concentrated in or near the cathode and when they have not had sufficient time to flow back. This is indeed what was observed for actuation at lower electrical voltages where the measurements were taken after 10 min of actuation (for SD) and 10 min of short-circuit (for RSD) and the operation was done for the next voltage right after. This does not allow sufficient time for the IL to flow back (10 min of short-circuit compared to 20 h) and results every time in a significant RSD varying between 0.6% and 0.9% (Figure 4c). The energy supplied then has the simple role of reswitching part of the cations that had flowed back during the period when the actuator was short-circuited. In the configuration where the cations remain close to the cathode, then the SD is proportional to the electrical voltage supplied. Figure 4b also shows that if the experiment is repeated under the same conditions after 20 h of rest, the SD was repeatable provided that the RSD before each electrical voltage applied was taken into account.

The RSDs measured after the short-circuits as a function of the various applied voltages are shown in Figure 4c. For both sets of measurements, a decrease in RSD was observed as the applied voltage decreased. This indicates that the actuator returns more easily to the vertical position as more electrical voltages are applied one after another. This was also confirmed by the lower RSDs of the second series compared to the first. The decrease in RSD with the successive application of electrical voltages, which indicates that the beam returns to a position increasingly closer to vertical after a short circuit, could be explained by an increase in the mobility of the ions which is facilitated the more electrical voltages are applied.

CONCLUSION

In this paper, the strain and its dynamics of PEDOT:PSS polymer-based trilayer micro-actuators were studied. The measurements were carried out by applying a DC voltage for 10 min followed by a short circuit for 10 min before any other operation. The polarity was preserved between each measurement; only the amplitude of the signal was modified. The SD measured after applying a voltage for 10 min and the RSD measured after a 10 min short circuit were thoroughly analyzed. This study has highlighted the following points: (1) The purely mechanical relaxation of a micro-actuator is much faster than that of an actuator previously subjected to an electrical voltage and a short circuit; (2) In the latter situation the actuator returns to its vertical position after several days, but after 10 min an RSD could be measured for a relatively stable state of the micro-actuator; (3) The RSD observed indicates that the cations did not return to their original position after the short-circuit; and that after 10 min, the reflux of the ions was slow; (4) The overall SD obtained depends on the previous RSD. Then, the position of the beam in space, when it is energized, or after a short circuit, varies slightly depending on the previous states imposed on this beam; (5) The first electrical voltage applied to the actuator is more efficient than the subsequent ones, favored by a prior homogeneous cation distribution in the material; (6) Without time in standby, the SDs obtained successively evolve linearly with the voltages applied as long as the previous state of the beam (RSD) is taken into account; (7) The successive application of a DC electrical voltage with the same polarity promotes an increase in the mobility of the ions, characterized by a continuously decreasing RSD. These phenomena, which are directly related to the intrinsic operation of the actuator, are similar to a memory effect resulting from the combination of the distribution of the ions in the actuator and their low mobility.

Figure 1 :

 1 Figure 1: Schematic representation of the volume variation mechanism induced by the movement of ions in PEDOT:PSSbased trilayer actuators.

Figure 2 :

 2 Figure 2: a) Chronovoltamogram of the measurement; b) PEDOT:PSS-based actuator in its initial state, before the first electrical connection; c) After application of a voltage of 1.75 V for 10 min; d) Picture taken of PEDOT:PSS electrodes after a short-circuit of 10 min to measure a residual strain difference (RSD1); e) Second application of a voltage of 1.75 V for 10min; f) Picture taken after a 10 min short-circuit to measure the residual strain difference (RSD2).

Figure 3 :

 3 Figure 3: a) Chrono-voltamogram of a voltage supply followed by a short-circuit, and b) normalized (divided by the maximum value) RSD determined from the curvature radius as function of time after an applied voltage of 1.75 V for 10 min followed by short-circuit. The relaxation was recorded for 25 h. c) RSD of the purely mechanical relaxation of a beam containing IL recorded for 35 min: manual strain maintained for 10 min and then released.

Figure 4 :

 4 Figure 4: a) Chronovoltamogram of the two series of SD measurements. b) SD as function of decreasing applied voltage. c) RSD of the microbeam after a short-circuit as function of the applied voltage.

ACKNOWLEDGMENTS

This work was supported by the H2020 project TWINNIMS (Grant agreement 857263), the French Government through the National Research Agency (ANR) and the MicroTIP, ROBOCOP, PIA EQUIPEX LEAF (ANR-11-EQPX-0025) projects, and the French RENATECH network.


N Hadda

L Sibille

F Nicot

R Wan

F Darve

A microstructural cluster-based description of diffuse and localized failures

This paper presents the analysis of microstructural mechanisms observed during both localized and diffuse failures in granular media, highlighting similarities and differences in their characteristics. Twodimensional DEM (Discrete Element Method) granular assemblies with medium dense and dense packings were subjected to different biaxial loading paths to induce either a localized or diffuse failure mode. A cluster based analysis is proposed to investigate interactions of particles at the mesoscale through the calculation of the second-order work from microscopic variables. It is shown that such analysis of clusters together with second order work facilitates precise shear band pattern recognition during localized failure. The evolution of such clusters in terms of their spatial distribution, size and number of particles involved, as well as the role played by strong and weak phases, describes the nature of the failure at every stage during loading. Such microstructural descriptors can predict the propensity of the specimen to fail either according to a diffuse or a localized mode.

INTRODUCTION

The identification of the underlying micromechanics of failure processes in granular media remains elusive despite numerous studies on the various factors leading to failure. Failure has always been related to either an instability criterion such as the vanishing of the second order work [START_REF] Hill | A general theory of uniqueness and stability in elastic-plastic solids[END_REF] or associated with contact networks (Radjai et al. 1999, Iwashita and[START_REF] Iwashita | Micro-deformation mechanism of shear banding process based on modified distinct element method[END_REF] namely strong and weak phases, connectivity, cycles and force chain networks [START_REF] Tordesillas | Evolution of functional connectivity in contact and force chain networks: Feature vectors, k-cores and minimal cycles[END_REF], force chain buckling [START_REF] Tordesillas | Force chain buckling, unjamming transitions and shear banding in dense granular assemblies[END_REF]) and inter-particles relative motions [START_REF] Kuhn | Contact rolling and deformation in granular media[END_REF].

Recently, it was shown [START_REF] Nicot | Inertia effects as a possible missing link between micro and macro second-order work in granular media[END_REF]) that the genesis of failure in an assembly of particles is closely related to the sign of the second order work at contacts and its spatial evolution during loading history [START_REF] Hadda | Micromechanical analysis of second order work in granular media[END_REF]. For instance, based on directional analyses [START_REF] Gudehus | A comparison of some constitutive laws under radially symmetric loading and unloading[END_REF]) performed on a dense and loose DEM specimens, it was found that the macroscopic second order work is in good agreement with the second order work computed from microscopic variables. Furthermore, the same study revealed that local instabilities within a granular assembly described by negative or zero values of the second order work at the contact level are in line with the vanishing of the macroscopic second order work. This finding was obtained from some statistical analyses focusing on the population and the spatial distribution of contacts exhibiting negative values of second order work, herein denoted c -, along all strain probe directions, including those which were proven to be unstable.

For instance, the spatial distribution of c -contacts is closely related to the emergence of the failure mode, i.e. diffuse and localized. As such, it seems that the understanding of failure process and the triggering of modes is closely associated with how disperse or concentrated the spatial distribution of c -contacts is. The concentration of c -contacts conveys the notion of proximity between unstable particles and in the limit refers to finding the most loaded particles with the most c -contacts. Hence the notion of local instability is most suitably pursued by calculating the sum of second order work at all c -contacts for a given particle. In this paper, and in continuity with our previous studies, the agreement between microscopic and macroscopic second order works is briefly reviewed along a drained biaxial compression. Then, attention is focused on the extent to which c -contacts relate to the observed failure mode, namely localized and diffuse. Finally, a cluster based approach is introduced in order to highlight the key role of c -contacts in describing failure at the mesoscale through an aggregation phenomenon.

FAILURE MODE ALONG DIFFERENT LOADING PATHS

Numerical model

It is well recognized that the failure mode depends essentially on the density of the granular packing and may be influenced by the loading direction. Thus, in order to achieve distinct failure modes through numerical simulations, two specimens S 1 (dense) and S 2 (medium dense) were considered and each was subjected to a different loading in order to undergo a localized and a diffuse failure respectively. Both specimens S 1 and S 2 originate from the same two-dimensional discrete element model [START_REF] Cundall | A discrete numerical model for granular assemblies[END_REF]) consisting of 21,000 particles and they differ only in their initial density (packing). The interparticle interaction is governed by a cohesionless contact law and involves three mechanical parameters: the normal contact stiffness k n , the tangential contact stiffness k t and the friction coefficient µ incorporated at the contact level through the Coulomb friction law. The characteristics and mechanical parameters are given in Table 1.

The open source discrete element code Yade (Šmilauer et al. 2010) was used to perform all the simulations presented in this paper.

The granular assemblies were both compacted from an initially sparse cloud of particles to an isotropic state by moving the four rigid frictionless walls toward the center of the specimen until the isotropic pressure Figure 1. Evolution of the deviatoric stress q and the volumetric strain ε v in terms of the axial strain ε 1 along the drained biaxial loading performed on S 1 and S 2 .

σ 1 = σ 2 = 300 kPa is reached. Both specimens S 1 and S 2 were then subjected to a drained biaxial compression. This compression is carried out by imposing a constant strain rate (ε 1 = 0.01 s -1 ) as the walls move in the axial direction, and maintaining constant the lateral pressure using a servo-control. In a second simulation, the medium dense specimen S 2 was loaded along a proportional strain path. This type of loading constraints the lateral strain to evolve proportionally to the axial strain at each time step such that ε 2 =-ε 1 /R (where R is a constant) as walls move in the axial direction at a constant strain rate (ε 1 = 0.01 s -1 ). Note that for R < 1, the loading path is dilatant, for R > 1, the loading path is contractant and the case of R = 1 coincides with the isochoric loading path. In this study, S 2 was loaded along a dilatant proportional strain path with R = 0.8.

The evolutions of the deviatoric stress q = σ 1σ 2 and the volumetric strain ε v in terms of the axial strain ε 1 for S 1 and S 2 resulting from the classic drained compression are plotted in Figure 1. As described, the specimen S 1 shows a typical behavior of a dilatant dense granular medium with a distinctive q peak, while the specimen S 2 exhibits a behavior akin to a medium dense granular medium with a volumetric contractancy up to 5% of axial strain followed by dilatancy.

Localized failure

According to the incremental deviatoric strain ε dev fields displayed in Figure 2-a, the dense specimen has clearly underwent a localized failure with a clear shear band after the limit stress state is reached. Incremental volumetric strains, ε dev , were computed between each successive pair of stress states identified in Figure 1. At the beginning (see snapshot 1), the strain field is Figure 3. Evolution of the deviatoric stress q and σ 1σ 2 /R in terms of the mean pressure p and the axial strain ε 1 respectively along the proportional strain path loading performed on S 2 .

homogenous. The snapshot 2 shows not only the fact that strain localizations are well developed at the peak, but also these have already started to define the shear band pattern observed after peak through several thin incipient shear bands. After the peak (see snapshots 3-5), a unique shear band is observed while all others have disappeared. Contrary to what has been observed for S 1 , it appears that there is no persistent shear band appearing after failure for the medium dense specimen S 2 as shown in Figure 2-b. However, the incremental deviatoric strain fields computed before peak appear to be similar for both S 1 and S 2 samples.

Diffuse failure

Based on these results reported in the previous section, the medium dense specimen, compared to the dense one, is more susceptible to undergo a diffuse failure along a proportional strain path with a reasonable value of R. Figure 3 shows the evolution of the deviatoric stress q in terms of the mean pressure p. The deviatoric stress q reaches a peak and then both q and p start to decrease until they both vanish. The peak of q for such loading path is not of relevancy when it comes to failure; rather one should consider the peak of σ 1σ 2 /R which represent the mixed limit stress state at which the second order work criterion is verified [START_REF] Nicot | Granular media failure along triaxial proportional strain paths[END_REF]. Figure 3 shows that σ 1σ 2 /R also passes through a peak, which is slightly exceeded afterward before both σ 1 and σ 2 vanish.

The incremental deviatoric strain fields corresponding to the stress states chosen along the proportional strain path displayed in Figure 4 seem to be quite well homogenous until the vanishing of stresses. Very weak strain localizations can show up intermittently, but they soon fade away. In this case, and on the contrary of what has been observed along the drained compression, the failure is diffuse and develops within the whole granular packing instead of being localized within a shear band.

SECOND ORDER WORK FROM MICROSCOPIC VARIABLES

For granular media, it has been shown [START_REF] Nicot | Inertia effects as a possible missing link between micro and macro second-order work in granular media[END_REF]) that the second order work W 2 can be expressed in terms of microscopic variables as:

The first term of Eq. ( 1) is a summation over contacts and involves the incremental variations of the contact force f c and branch vector lc . The second term, a summation over particles, accounts for the inertial effects that may occur during loading. It involves the incremental variation of the resultant force f p applied on the particle and the variation in position

x experienced by the particle during loading.

In the absence of inertial effects (quasi-static regime), which is the case in this study, the contribution of the second term is negligible and W 2 is reduced to:

First of all, it is interesting to compare the second order work computed from microscopic variable W m variable W 2 = σ : ε [START_REF] Nicot | Bifurcation and second-order work in geomaterials[END_REF], where σ and ε are the Cauchy stress tensor and strain tensor respectively, and to see to what extent they are equal along a drained biaxial compression. For this purpose, both second order works were computed between each two successive stress-states chosen at a regular interval along the loading path.

Figure 5 shows that the two second order work expressions, computed differently, agree with each other especially before failure along the classic biaxial compression performed on S 1 . The vanishing of the second order work corresponds to the peak of q (i.e. the limit stress state). After failure, a deviation between the two curves may be observed, especially when the specimen undergoes a sharp softening (local avalanche). This phenomenon is usually accompanied by the release of inertial effects, whose contribution may become effective, and thus should be accounted for while computing the second order work from microscopic variables.

This good agreement is confirmed through the spatial distributions of c -contacts identified at the stress-states 1-5 (Figure 6-a) and 1-4 (Figure 6-b) along the biaxial compression performed on S 1 and the proportional strain path performed on S 2 respectively. Let's recall here that c -contacts account for those which exhibit a zero or negative value of w 2 = f c lc .

Each c -contact is presented by a black dot in the snapshots displayed in Figure 6. At the beginning, the spatial distribution of c -contacts along the biaxial compression is quite homogenous (snapshots 1 and 2). Further c -contacts become more and more concentrated within the shear band (characterized by the incremental deviatoric strain field) as the axial strain increases. However, the width of the shear band defined by c -contacts is relatively larger than the one derived from the incremental deviatoric strain field. For the proportional strain loading path, c -contacts seems to be homogeneously distributed before and after failure with no particular distinguishable structure. Note that in the last snapshot, the sparse number of c -contacts is due to the loss of contacts within the sample resulting from the static liquefaction, and not to a switch of contacts from negative to positive values of w 2 .

The second order work computed from microscopic variables was proven to be in line with the vanishing of the macroscopic second order work. Moreover, the mode of failure can be directly predicted by means of the spatial distribution of c -contacts. Hence, it definitely constitutes a powerful and suitable tool for the identification of the underlying micro structural mechanisms governing the failure process.

CLUSTERS DEFINITION BASED ON c-CONTACTS

Based on the results shown above regarding c -contacts evolution and distribution during diffuse and localized failures, and according to other findings [START_REF] Hadda | Micromechanical analysis of second order work in granular media[END_REF] reporting the strong relation between the increase of c -contacts number and instability, the distance between c -contacts seems to play a key role in triggering failure within a granular assembly, especially when the latter come close together. It is definitely clear, that c -contacts concentrated within the shear band (Figure 6-a) are much closer to each other compared to those outside the band. Moreover, during diffuse failure, the proportional strain loading path seems to prevent c -contacts from getting much closer, dispersing their concentrations before they get larger. Thus, it is quite interesting to focus the investigations on c -contacts gatherings for a better understanding of failure triggering and process. When c -contacts come so close to each other, this results into another or more contacts belonging to a given particle turning into c -contacts. Thus, a gathering of c -contacts can be identified from the most loaded particles and such particles, are most likely to imply as they sustain more and more c -contacts.

Regardless of the number of c -contacts compared to the total number of contacts belonging to a given particle, if the condition (3) is fulfilled, then the particle is likely to be part of a cluster made of at least two particles in contact of that kind. Hence, such particle is denoted by p -.

In this section, interest will be focused on the evolution and distribution of p -clusters along both considered loading paths which led to different failure modes.

Figure 7 shows respectively the evolution of the total number, the mean size of p -clusters and the percentage of particles involved within the clusters for S 1 while undergoing a localized failure. The largest number of clusters formed is reached at the peak of q, involving naturally the largest number of particles. Then, a sharp decrease is observed just after the peak followed by a fluctuation around a constant value.

The ongoing increase of the mean size of p -clusters during the sharp decrease of the number of clusters right after the peak emphasizes a crucial phase of the localized failure process. This phase corresponds to the appearance of the persistent shear band, throughout which, already existing clusters enlargement and extinction processes take over cluster growth. This result is consistent with what has been observed through the evolution of the incremental deviatoric strain during localized failure (Figure 2-a).

Identically, Figure 8 shows p -clusters evolution for S 2 while undergoing a diffuse failure. Surprisingly, the evolution of the number of p -clusters does not show any particular tendency of failure occurrence. In contrast with localized failure, diffuse failure is characterized by an almost constant number of p -clusters all along the proportional loading path after a rapid modest increase at the beginning of the loading and a spontaneous extinction when static liquefaction fully took place.

Figure 8 shows also a very narrow spread of the mean size of p -clusters compared to the one observed during the localized failure, which means that p -clusters formed during such loading preserve almost the same size with a very restricted extension. However, at each peak of σ 1σ 2 /R, a slight increase of p -clusters mean size is observed with almost the same number of particles involved within clusters before the peak. Thus, some clusters grow larger during these failure stages but obviously their enlargement is very limited in size and very short in time.

Noticeably, the total number of p -clusters formed during diffuse failure overwhelm those formed during localized failure and involve more p -particles. More or less 30% of particles are involved within p -clusters during diffuse failure while this value is reached only at peak during localized failure and is further reduced to 15%. This phenomenon, which persists during diffuse failure, occurs only at peak during localized failure.

According to Figure 9, where spatial distributions of p -clusters are displayed at stress states chosen along both loading paths, p -cluster distribution proves to be more suitable than c -contacts to directly describe the observed failure mode. For instance, for localized failure, large clusters are found to be distributed within shear bands. Also, the width of the shear band characterized by p -clusters is not as wide as the one obtained with c -contact distribution (Figure 6-a), but it rather recalls the one observed in the incremental strain fields (Figure 2-a). Moreover, the p -cluster spatial distribution was able to point out the slight strain localizations occurring during diffuse failure, which are hardly distinguishable through c -contact distribution.

CONCLUSIONS

The present study provides a method by which the origins of localized and diffuse failure modes can be deciphered through the analysis of particle clusters in combination with second order work computations. As pointed out in the paper, the DEM studies on both medium dense and dense granular assemblies subjected to biaxial loadings paths lend considerable support to the idea that localization patterns and their evolution are driven inherently by the clustering of particles with contacts violating the second order work. More precisely, the consideration of so-called p -cluster particles for which the net second order work of all contacting particles is zero or negative leads to a richer interpretation of failure modes, i.e. either localized or diffuse.

Figure 2 .

 2 Figure 2. Incremental deviatoric strain fields computed at the stress states 1 to 5 along the drained biaxial compression performed on S 1 (a) and S 2 (b).

Figure 4 .

 4 Figure 4. Incremental deviatoric strain fields computed at the stress states 1 to 4 along the proportional strain path loading performed on S 2 .

Figure 5 .

 5 Figure 5. Comparison between the second order works computed from both macroscopic and microscopic variables along the drained compression and the proportional strain loading path performed on S 1 and S 2 respectively.

Figure 6 .

 6 Figure 6. Spatial distributions of c -contacts computed at the stress states 1 to 5 along the drained compression performed on S 1 (a) and states 1 to 4 along the proportional strain loading path performed on S 2 (b).

Figure 7 .

 7 Figure 7. Evolution of the number, mean size of p -clusters and the number of particles involved within clusters along the drained compression performed on S 1 .

Figure 8 .

 8 Figure 8. Evolution of the number, mean size of p -clusters and the number of particles involved within clusters along the proportional loading path performed on S 2 .

Figure 9 .

 9 Figure 9. Spatial distribution of p -clusters at the stress states 1 to 5 along the drained compression performed on S 1 (a) and states 1 to 4 along the proportional strain path loading performed on S 2 (b).

Table 1 .

 1 Physical and mechanical parameters of both specimens S 1 (dense) and S 2 (medium dense).

		S 1	S 2	Unit
	Particles diameter range Density (ρ)	6-18 3,000		mm kg/m 3
	k n /D(*) k t /k n µ	356 0.42 0.7		(MPa) --
	κ = kn p (**) (p = 300 kPa)	1,200		-
	Height / Width	1.85	1.35	-
	Void ratio (e)	0.174	0.217	-
	Coordination number (z)	4.20	3.36	-

*D denotes the mean diameter between two particles in contact. **κ and p denote respectively the normalized contact stiffness and the mean pressure.

with the second order work computed from tensorial


Alexander Mendler
email: alexander.mendler@tum.de

Michael Döhler
email: michael.doehler@inria.fr

Christian U Grosse
email: grosse@tum.de

SELECTION OF DAMAGE-SENSITIVE FEATURES BASED ON PROBABILITY OF DETECTION CURVES

Keywords: Structural health monitoring, ambient vibrations, Fisher information, probability of detection

The first phase of each structural monitoring project is the operational evaluation. Its purpose is to define relevant damage mechanisms and environmental conditions, to consider the data acquisition limitations on site, and to justify the investment. Subsequently, relevant measurement quantities and damagesensitive features are selected, but very few systematic approaches exist in the literature on how to select the most appropriate features. The presented paper fills this gap and develops an approach to select damage-sensitive features based on probability of detection (POD) curves. The POD curves are generated based on a novel method for statistical damage detection tests that requires a finite element model and vibration data from the undamaged structure. However, no data is required from the damaged state, making it particularly suited for unique or large and complex engineering structures. The approach explicitly considers the uncertainties in the features due to unknown loads, measurement noise, and short measurement durations. Although global damage-sensitive features are considered, such as modal parameters and subspace-based residuals, the detectability is evaluated for local structural components. The paper includes a proof of concept study on a laboratory structure. The results demonstrate that the developed method successfully finds the feature with the highest damage detectability for a chosen damage scenario, and that the detectability varies depending on the monitored local component.

INTRODUCTION

National building standards require visual inspections and non-destructive testing for critical engineering structures in periodic intervals. Increasingly, structural health monitoring (SHM) system are installed, meaning sensors are permanently installed on the structure and online damage diagnosis algorithms are trained to automatically diagnose damages in-between the scheduled inspections. The four main stages of each SHM project are the operational evaluation, the data acquisition, the extraction of damagesensitive features (such as natural frequencies or mode shapes), and their statistical evaluation [START_REF] Farrar | Structural health monitoring: A machine learning perspective[END_REF]. Each phase exhibits distinct challenges, for example, the correct parametrization of damage, the selection of appropriate measurement quantities and sensor locations, and the selection of appropriate damagesensitive features, and algorithms.

Feature selection is a critical topic. Some machine learning experts prefer a 'blind' feature selection, for example, based on neuronal networks whose training process can involve the automated selection of features. Others prefer a 'manual' selection based on engineering judgement, experiences with similar structures, or analytical and structure-specific approaches. Such analytical approaches include sensitivity analyses using finite element models [START_REF] Farrar | Structural health monitoring: A machine learning perspective[END_REF], where the effect of individual parameter changes on the data-driven features is studied, or approaches based on confidence intervals, where the changes in the features is put into relation to the statistical uncertainties that are inherent to the features' estimation process. However, few method-specific approaches exist to select the feature with the highest damage detectability. A novel idea is to select features based on probability of detection (POD) curves for critical damage scenarios, as they quantify the reliability of a SHM system. POD curves are standard approaches to verify the effectiveness of non-destructive testing [START_REF] Berens | NDE Reliability Data Analysis-Metals Handbook[END_REF], but they are rarely used in SHM. The reason for this is that most methods (e.g., the 29-29 method, the â vs. a method, the hit/miss method, the delta method, Bayesian approaches) require empirical data from the damaged state, and this data is typically not available for large and unique structures, such as bridges, highrises, and dams. The main idea of this article is to propose a new method for the creation of POD curves that does not require empirical data from damaged structures. The second objective is to apply this approach to various features and to select the feature that exhibits the highest damage detectability.

The paper is organized as follows: Section 2. recaps state-of-the-art vibration models. Section 3. explains how damage-sensitive features can be formed based on the dynamic system properties, how they can be evaluated statistically, and how POD curves can be predicted based on measurement data from undamaged structures. Section 4. showcases a proof of concept study based on a laboratory beam with extra masses, and Section 5. summarizes the main findings.

BACKGROUND

The presented approach is applied to vibration-based features and requires in-depth knowledge on operational modal analysis. This section recaps how modal parameters can be obtained from measurement data using subspace-based system identification. All following considerations are based on linear and time-invariant dynamic systems with m degrees of freedom (DOF)

Mü + D u + Ku = f, (1)
where M, D and K are the mass, damping, and stiffness matrices R m×m , u ∈ R m is the displacement vector, and f ∈ R m is the force vector. Sampling the displacement vector at ∆t and substituting the state vector x k = u(k∆t) T u(k∆t) T T ∈ R n , the model is transformed to a discrete-time state space model, with n = 2m,

x k+1 = Ax k + w k y k = Cx k + v k (2)
where A ∈ R n×n is the state transition matrix, C ∈ R N ch ×n is the output matrix and N ch is the number of measurement channels. The term y k is the output vector, and w k and v k state noise and measurement noise.

[3]

Subspace Decomposition

The first step is to acquire output data Y = [y 1 y 2 . . . y N ] T ∈ R N ×N ch (displacements, velocities, or accelerations), where N is the number of data points and N ch is the number of measurement channels. Next, output covariance functions are evaluated

R i = 1 N -i N -i
k=1 y k+i y T k and arranged in a block Hankel matrix, i.e., a matrix with identical blocks on the anti-diagonals

H p+1,q =      R 1 R 2 . . . R q R 2 R 3 . . . R q+1 . . . . . . . . . . . . R p+1 R p+2 . . . R p+q      = U 1 U 0 S 1 0 0 S 0 V T 1 V T 0 , (3)
where p and q are time lag parameters. Thirdly, singular value decomposition (SVD) is applied to the block Hankel matrix, and the resulting quantities are truncated at the presumed model order, determined through the m physical modes of vibration, yielding the first n = 2m singular values in S 1 and the corresponding singular vectors U 1 (the column space). The remaining singular values S 0 contain noise with the corresponding left null space vectors in U 0 .

Modal Identification

This section revisits how modal parameters can be derived from the block Hankel matrix [START_REF] Van Overschee | Subspace identification for linear systems: Theory, Implementation, Application[END_REF]. First, the observability matrix O p+1 is constructed based on the first n singular values and the column space

O p+1 = U 1 S 1/2 1 . (4)
The output matrix C can be extracted from the observability in Eq. ( 4) as the first block row, and the state transition matrix can be approximated through regression, using the shift-invariance property [4]

A = (O ↑ p+1 ) † O ↓ p+1 , (5)
where O ↑ p+1 and O ↓ p+1 denote the observability matrix without the last and without the first block row, respectively. Modal parameters can be estimated by solving the eigenvalue problem AΦ = ΦΛ. Eigenvalues occur in complex conjugate pairs Λ = Φ -1 AΦ = diag(λ 1 , λ1 . . . , λ n , λn ) and so do the eigenvectors. Ultimately, natural frequencies and mode shapes are calculated as follows

µ j = log λ j ∆t , f j = |µ j | 2π , Ψ = CΦ = Ψ 1 • • • Ψ n . (6)

Mode Shape Normalization

Mode shapes are unique dynamic properties up to a scaling factor. To enable a comparison of mode shapes, they have to be normalized, e.g., to unit displacement of one of its components l,

Ψi = 1 max{|Ψ i,l |} • Ψ i . (7)
This forces the maximum amplitude of each normalized mode shape to be equal to one (even when damage has occurred). Consequently, the maximum mode shape component carries no damage-related information after scaling and the mode shape dimension might as well be reduced to Ψred ∈ R N ch -1×m , which is done in this paper.

METHODOLOGY

Data-driven Residuals

This section shows a selection of damage-sensitive residuals that can be formed based on dynamic system properties from the previous section. In the subsequent sections, the POD curves will be constructed for all features presented in this section. For example, the natural frequencies could be extracted and compared to the baseline values from the training phase

r 1 = f -f 0 , (8)
where the number of entries in the residual vector depends on the number of observed modes of vibration.

Another example is the mode shape residual. To reduce the mode shape matrix into a single vector, the vectorization operator vec(•) is used and all matrix columns are stacked

r 2 = vec( Φred ) -vec( Φred,0 ), (9)
where the number of residual entries now depends on the number of modes of vibration and the number of sensors. Since modal frequencies and mode shapes are the results from the same estimation procedure, they can (and should) also be analyzed simultaneously, using a modal vector η with the corresponding residual

r 3 = η -η 0 , where η = f vec( Φred ) . ( 10
)
The last considered residual is the subspace-based residual. It takes advantage of the orthogonality between the left null space singular vectors and the block Hankel matrix in the reference state [START_REF] Basseville | Subspace-based fault detection algorithms for vibration monitoring[END_REF]

r 4 = vec(U T 0 H p+1,q ) (11)
and its size depends on the time lag parameters p and q.

Damage Detection

This section summarizes the specific damage diagnosis method for which the POD curves are constructed. In this paper, structural changes are analyzed based on a parametrized statistical test, i.e., a statistical test on data-driven residuals that is linked to structural models through sensitivity vectors. The test statistic is defined as [START_REF] Basseville | Subspace-based fault detection algorithms for vibration monitoring[END_REF]

t = N • r T Σ -1 J J T Σ -1 J -1 J T Σ -1 r, ( 12
)
where N is the sample size during testing, r is one of the damage-sensitive residual vectors from Section 3.1., J is the first-order derivative of the residual with respect to structural design parameters θ

J = ∂E θ [r] ∂θ θ=θ 0 , ( 13
)
and Σ is the sample covariance of the damage-sensitive residual which, when evaluated on n b datasets of length N , is defined as

Σ = N n b -1 n b k=1 r k r T k . ( 14
)
Structural parameters θ ∈ R H typically include material properties, cross-sectional values, prestressing forces, or geometric properties of the structure that manifest damage. Defining the relevant damage parameters is part of the first phase of each monitoring phase, the operational evaluation.

The Jacobian matrix J links structural parameter changes and changes in data-driven features through r ≈ J θ. It can either be approximated using the finite difference method, or it can be calculated analytically. For the subspace-based residual, the analytical calculation is described in [START_REF] Allahdadian | Towards robust statistical damage localization via model-based sensitivity clustering[END_REF] and for the modal residuals, typical approaches include the modal approach [START_REF] Fox | Rates of change of eigenvalues and eigenvectors[END_REF] or Nelson's method [START_REF] Nelson | Simplified calculation of eigenvector derivatives[END_REF]. The covariance matrix is calculated based on data, where the main diagonal holds the variance of each residual entry. It can either be computed empirically based on multiple sets of features from the training phase, or in a combined analytical/empirical approach based on a single measurement record for both the modal residuals [START_REF] Reynders | Uncertainty bounds on modal parameters obtained from stochastic subspace identification[END_REF][START_REF] Döhler | Efficient multi-order uncertainty computation for stochastic subspace identification[END_REF] and the subspace-based residuals [START_REF] Basseville | Subspace-based fault detection algorithms for vibration monitoring[END_REF][START_REF] Viefhues | Statistical subspace-based damage detection with estimated reference[END_REF].

Uncertainties due to unknown loads and measurement noise cause the test statistic in Eq. ( 12) to be scattered. Since the damage-sensitive residuals can be approximated by a Gaussian distribution, the tests can be approximated as a chi-squared distribution χ 2 (ν, λ) with ν = rank(J T Σ -1 J ) degrees of freedom, and non-centrality λ in the damaged state. To make a decision on whether or not the structure is damaged, a safety threshold t crit is introduced, for example, based on the distribution in the undamaged state, and an alarm is issued if the test exceeds the safety threshold, see Fig. 1.

Probability of Detection Curves

A powerful aspect of the parametric hypothesis tests in Eq. ( 12) is that the mean test response to a welldefined parameter change can be predicted based on quantities that are available in the undamaged state of the structure. If damage is assumed to be restricted to a single parameter change θ h -θ 0 h , the mean test response (the non-centrality) λ is [START_REF] Mendler | A reliability-based approach to determine the minimum detectable damage for statistical damage detection[END_REF]

λ = N (θ h -θ 0 h ) 2 F hh , F hh = J T h Σ -1 J h , (15)
where F hh is the Fisher information and J h is the column in the Jacobian matrix that corresponds to the examined parameter θ h . The larger the mean test response from Eq. ( 15), the more frequently the test statistic yields values beyond the safety threshold, and the higher the POD. The POD can be quantified as the area under the probability density function of the test statistics beyond the safety threshold, Fig. 1,

POD = ∞ t crit f χ 2 (ν, λ)(t)dt, (16)
where the safety threshold could, for example, be defined based on the quantile value of the healthy state distribution of the test statistic. Assuming that the test response λ is known, the formula from Eq. (15) can also be solved for the minimum detectable parameter change

∆ h = 1 θ h λ N F hh . (17)
Utlimately, the POD curves can be constructed as follows: First, the non-centrality is fixed to zero λ = 0.

Next, the non-centrality is gradually increased λ = λ + ε in steps of ε, while evaluating the POD based on the corresponding probability density function in Eq. ( 16) and Fig. 1 and the minimum detectable damage ∆ h using Eq. ( 17). Drawing the POD over ∆ h yields the probability of detection curve for parameter θ h .

APPLICATION

The goal of this section is to analytically construct POD curves based on data from an undamaged structure. Subsequently, actual damage is applied, and it is demonstrated that the empirical POD is close to the analytical one, which was predicted based on undamaged data.

Laboratory Experiment

The specimen under consideration is a hollow structural steel beam (HSS 152x51x4.78 mm) with a modulus of elasticity of E = 210,000 MPa and a total mass of m = 56.8 kg, Fig. 3. The beam is located at the University of British Columbia and has already been described in the referenced literature [START_REF] Mendler | A reliability-based approach to determine the minimum detectable damage for statistical damage detection[END_REF][START_REF] Mendler | Sensor placement with optimal damage detectability for statistical damage detection[END_REF].

The length of the beam is 4.1 m and pin supports are installed on both ends. The instrumentation consists of eight vertical accelerometers with a weight of 1.28 kg, placed at equal distances, and one shaker with a mass of 3.6 kg and a moving mass of 360 g that injects white noise signals in-between sensors 2 and 3. All other signal processing parameters are given in Table 1.

Analytical Results

All studies in this section are based on measurement data from the undamaged structure, where only data from the first four measurement channels is used, Fig. Sensor location The Fisher information F can be calculated based on the Jacobian matrix J and the covariance matrix Σ, see Eq. ( 15). The covariance matrix is calculated based on [START_REF] Basseville | Subspace-based fault detection algorithms for vibration monitoring[END_REF][START_REF] Döhler | Efficient multi-order uncertainty computation for stochastic subspace identification[END_REF] and for visualization, the standard deviations of the measured modal parameters (i.e., the square root of the covariance's main diagonal values) are shown in Fig. 4 through error bars. For modal parameters, the Jacobian matrix computation can be done based on the numerical model but for the subspace-based residual, data-driven components have to be considered [START_REF] Allahdadian | Towards robust statistical damage localization via model-based sensitivity clustering[END_REF]. Once the Jacobian and the covariance matrix are estimated, the Fisher information can be calculated and the main diagonal values F hh can be extracted.

The last required quantity is the safety threshold value t crit . It is different for each damage-sensitive residual, because the χ 2 -distribution exhibits a different number of degrees of freedom. In this study, the thresholds are set based on the 0.1% quantile value of the healthy state distribution, meaning one out of 1,000 test statistics is beyond the safety threshold in the undamaged state. The corresponding safety thresholds are 22 for the subspace-based residual, 16 for the frequency-based residual, 16 for the mode shape-based residual, and 22 for the combined frequency and mode shape-based residual.

With the Fisher information and the safety threshold at hand, the POD curves can be predicted for all damage-sensitive residuals from Eq. ( 8)- [START_REF] Viefhues | Statistical subspace-based damage detection with estimated reference[END_REF] with the results being displayed in Fig. 5. The POD curves are drawn for each monitoring parameter, i.e., for parameter 8 and 9 from Eq. ( 18). The first observation is that changes in parameter 8 lead to a much higher POD than in parameter 9, meaning the damage detectability close to the support is small. According to the POD curves for parameter 8, the frequencybased residual is the most sensitive to changes in the eighth mass parameter, as a parameter change of ∆θ 8 = 5% already leads to a POD close to 87%. The mode shape-based residual, on the other hand, is the least sensitive and a parameter change of 5% leads to a POD of close to 0%, and the subspacebased residual leads to a POD of about 20%. It should be noted that the POD curves are created for particular damage scenarios, captured through the parameter vector in Eq. ( 18). It appears that the mode shapes do not change for mass changes in beam segment eight and nine, but for other scenarios, they may significantly contribute to the damage detectability and should therefore be considered.

Empirical Validation

In this section, data from the damaged structure is analyzed for the first time to validate the predictions that were made based on data from the undamaged state. The main idea is to apply a 5% mass increase to the eighth beam segment, simply by placing two washers on top of the beam between sensors seven and eight (see Fig. 3, right side), and to verify whether the empirically evaluated POD is close to the predicted ones from Fig. 5. To evaluate the POD empirically, the test statistic is applied to 100 data sets from the damaged structure, with a measurement duration of 10 s for each segment, and the relative number of tests beyond the safety threshold is counted.

For the frequency-based residual, the empirically evaluated POD is 81.1%, Fig. 6, which is very close to the analytical value of 87%. Likewise, the mode shape-based residual leads to a POD of 0% and the subspace-based residual exhibits a POD of 19.4%, as predicted in Fig. 5. This concludes the validation study, as the analytical POD curves predict the empirical POD for all residuals with sufficient accuracy.

Discussion

Safety threshold. A peculiar observation appears to be that, in some cases, the mixed frequency and mode shape-based residual leads to a lower POD than the frequency-based residual although more information is available from the structure, see Fig. 5. The reason for this is that the frequency-based residual exhibits a narrower distribution in the undamaged state and a lower safety threshold value of 16 instead of 22. Therefore, an equivalent mean test response leads to a higher POD.

Sensor layout. The study in the previous section concluded that the frequency-based residual is more sensitive to damage than the subspace-based residual. It is important to understand that this result is only valid for the considered sensor layout, because both the sensitivity and the covariance of extracted features depend on the sensor location. To emphasize this, Fig. 7 shows how the POD curves change if the number of channels is increased from 4 to 8. Consequently, the number of mode shape coordinates In other words, the POD curves are powerful means for sensor placement optimization. A systematic approach to optimize the sensor layout can be found in [START_REF] Mendler | Sensor placement with optimal damage detectability for statistical damage detection[END_REF].

Measurement noise. One of the advantages of the developed approach is that the uncertainties in the feature due to measurement noise is explicitly considered in the covariance matrix, without having to measure the signal-to-noise ratio. To demonstrate this, vibration data is numerically generated and different noise levels are added to the vibrations. Fig. 7 shows the POD curves for two cases, the mode shape-based residual with a measurement noise level that corresponds to 5% and 10% of the output signal's variance. For the damage scenario ∆θ 8 = 5%, the POD increases from about 0% to 12% due to decreased noise level. For other residual formulations considered in this paper, the POD changed as well but the mode shape-based residual appears to be particularly sensitive to measurement noise.

Spatial information.

Since POD curves can be evaluated for different structural parameters, they carry spatial information. However, they should not be considered probability of localization (POL) curves, because the statistical test from Eq. ( 12) is not capable of localizing damage. More advanced statistical tests can identify the parameter that has changed and an approach to determine POL curves can be found in the referenced literature [START_REF] Mendler | On the probability of localizing damages based on mode shape changes[END_REF].

CONCLUSIONS

The main contribution of this paper is an approach to construct probability of detection curves, based on statistical tests and a finite element model. Advantages over other approaches are that no data from the damage state is required and the probability of detection can be analyzed based on data from undamaged structures. Moreover, the POD curves can be drawn for different structural parameters, so they include spatial information.

Secondly, it was shown that the developed method for POD curves is a powerful means for feature selection, as features with the highest damage detectability can be chosen before damage occurs. In this paper, the approach is applied to natural frequencies, mode shapes, combined frequency and mode shapes formulations, and the subspace-based residual, but many more features could be incorporated in the framework.

A laboratory case study demonstrates that the POD curves are accurate, even for real data and noisy measurement environments. The main findings are that the POD depends on the structural parameter that is monitored, the sensor layout, and the measurement noise levels. The effect of measurement noise appears to be particularly pronounced for mode shape testing, highlighting the importance of uncertainty quantification in structural health monitoring in general.

An interesting future research topic is the extension to additional features as well as a thorough investigation of measurement noise properties. If measurement noise is modelled accurately, POD curves could be constructed based on numerically generated data and findings can be transferred to real structures.

Figure 1 :

 1 Figure 1: Distribution of the test statistic in the undamaged and damaged state

Figure 2 :

 2 Figure 2: Schematic depiction of the laboratory beam indicating the eight beam segments

Figure 3 :

 3 Figure 3: Laboratory HSS beam with eight sensors and one shaker (left) and applied extra masses (right)

Figure 4 :

 4 Figure 4: Uncertainty quantification for natural frequencies and mode shapes

Figure 5 :

 5 Figure 5: Predicted POD curves based on data from the undamaged state for extra masses on beam segment 8 (left) and beam segment 9 (right) for different residual formulations from Section 3.1

Figure 6 :

 6 Figure 6: Damaged state: Comparing the empirical POD for different residual formulations from Section 3.1

Figure 7 :

 7 Figure 7: POD curve based on measurement data for an increased number of measurement channels (left), and based on simulated data with varying levels of measurement noise (right).

Table 1 :

 1 Data acquisition and signal processing parameters

				Data Acquisition			Data Segmenting
				Measured quantity	acceleration No. of training segments 100
				Sampling frequency	200 Hz	No. of testing segments	100
				Reference data length	110 min	Samples/segment	2,000
				Testing data length		110 min	Duration/segment	10 s
				Modal Parameters		Subspace-based Residual
				System order n		14	System order n	14
				Time lags p			20	Time lags p	11
		0					15
		-20		Power spectral density Stable mode with standard deviation	14.8
							14.6		1	2	3	4	5	6	7	8
		-40					14.4
	Power [dB]	-80 -60					13.8 14 14.2	Model order
		-100					13.6
							13.4
		-120
							13.2
		-140					13
		0	20	40	60	80	100
				Frequency [Hz]

2.

In this state, one long measurement with

ACKNOWLEDGMENTS

This research study is funded by dtec.bw -Digitalization and Technology Research Center of the Bundeswehr. The vibration data was kindly provided by Dr. Carlos E. Ventura from the University of British Columbia, Canada.