Faissal Merzaki
email: faissal.merzaki@insa-rennes.fr

Senior Member, IEEE Philippe Besni

Mohamed Himdi
email: mohamed.himdi@univ-rennes1.fr

Xavier Castel
email: xavier.castel@univ-rennes1.fr

Maëlle Sergolle
email: maelle.sergolle@naval-group.com

Thierry Levavasseur
email: t.levavasseur@seribase.fr

Patrick Caldamone
email: p.caldamone@seribase.fr

Philippe Besnier
email: philippe.besnier@insa-rennes.fr

A Compact Double-Sided FSS Absorbing Wall for Decoupling 5G Antenna Arrays

Keywords:

The fifth generation (5G) of mobile communications gives birth to the evolution of the base stations. They progressively include arrays of antennas in the sub-6 GHz frequency band worldwide. Moreover, 5G takes advantage of the beamforming capabilities to gain power efficiency and higher throughputs. Within this framework, achieving high decoupling performance between any pairs of antennas and compactness at the same time become a major challenge. The present paper describes an advanced solution based on the interposition of a thin and low-profile absorbing wall operating in the 3.4-3.6 GHz frequency band. The absorption is provided by a frequency selective surface whose elementary pattern size is first optimized. Then, the compact absorbing wall that associates back-to-back pieces of absorber is interposed between elements of an antenna array. Moreover, the fabrication of both antenna array and decoupling device involve lightweight and low-cost composite technology. The results show that such a compact decoupling device achieves excellent decoupling capabilities while keeping up the antenna array performance. Index Terms-Antenna decoupling, frequency selective surface (FSS)

coupling effects between communication systems operating at various frequency bandwidths. Moreover, 5G communications support multiple-input multiple output (MIMO) transmissions to handle beamforming, as a principle for increasing data rate transmission. Coupling issues within the involved array of antennas may jeopardize the performance of MIMO transmissions that is optimum if input/output signals are uncorrelated. This article focuses on this specific issue.

The literature describes numerous strategies to decouple neighboring antennas sharing, at the moment, the same printed circuit board. We may classify these strategies in three main categories.

First, the patterning of antennas common ground plane may prevent the flowing of currents from one antenna to another. This set of solutions is based on the so-called "defected ground structure" technique. For instance, such a solution may consist in etching a thin slot into the common ground plane of an antenna array [START_REF] Zakaria | Design of reconfigurable defected ground structure (DGS) for UWB application[END_REF]- [START_REF] Habashi | Mutual coupling reduction between very closely spaced patch antennas using low-profile folded splitring resonators (FSRRs)[END_REF]. This pattern can be seen as a resonator precluding the current from flowing from one antenna to another through the ground plane. Therefore, it decreases the induced voltage difference at neighboring antennas and reduces the coupling effect, accordingly.

Second, coupling between antennas may be due to surface wave propagation into the substrate. Therefore, another method consists in positioning "electronic band gap" (EBG) structure between antennas as introduced in [START_REF] Goussetis | Tailoring the AMC and EBG characteristics of periodic metallic arrays printed on grounded dielectric substrate[END_REF]. In this case, the EBG structure consists of a metamaterial surface (over the substrate) that prevents the propagation of surface waves through diffraction or absorption effects [START_REF] Farahani | Mutual coupling reduction in patch antenna arrays using a UC-EBG superstrate[END_REF]. However, this technique suffers from its size in case of strict compactness requirements for antenna arrays. Based on the tailoring of artificial magnetic conductors (AMCs) [START_REF] Goussetis | Tailoring the AMC and EBG characteristics of periodic metallic arrays printed on grounded dielectric substrate[END_REF] allowing reflection with a zero-phase shift, a double-sided piece of AMC surfaces was used as a vertical wall between monopole antennas [START_REF] Bayraktar | Composite planar doublesided AMC surfaces for MIMO applications[END_REF]. This solution applies to nonplanar antennas, but remains very useful for designing low-profile and narrow band antennas [START_REF] Feresidis | Artificial magnetic conductor surfaces and their application to low-profile highgain planar antennas[END_REF]. However, surface waves may have significant contributions for substrate with high dielectric constant. In the case of moderate dielectric constant, this is not a predominant effect.

Third, coupling also occurs at razing incidence through electromagnetic field propagation (in the near-field zone) over the air. A lot of research works have contributed to this with various concepts of thin absorbing materials including metamaterials [START_REF] Landy | Perfect metamaterial absorber[END_REF]- [START_REF] Kaur | Ultrathin wideband polarization independent compact metamaterial microwave absorber[END_REF] and frequency selective surface (FSS) structures [START_REF] Syed | A single-layer frequency selective surface for ultrawideband electromagnetic shielding[END_REF]- [START_REF] Costa | An equivalent circuit model of frequency selective surfaces embedded within dielectric layers[END_REF]. Beyond their intrinsic properties, their use as decoupling devices is illustrated in many applications [START_REF] Wang | A meta-surface antenna array decoupling (MAAD) method for mutual coupling reduction in a MIMO antenna system[END_REF]- [START_REF] Song | Suppression of the mutual coupling between two adjacent miniaturized antennas utilizing printed resonant circuits[END_REF]. These solutions mainly consist of implementing an absorbing structure as a wall between antennas to prevent the electromagnetic field propagation over the air and thereby to enhance the decoupling without significant modification of the antenna radiation pattern. Although these various techniques exist, they focus on decoupling high-profile antennas (such as monopoles). Only a few of them deal with the specific planar MIMO configurations for 5G applications and use suspended walls or superstrates. The solution investigated here is rather based on the design of a limited size decoupling wall for decoupling an array of patch antennas.

In this article, we present a low-cost, low-profile decoupling structure for a 5G planar antenna array, including three main new features. First, the decoupling structure is made of two back-to-back pieces of FSS absorbers forming a low-profile vertical wall with no significant impact on the antenna array performance. Second, it is fabricated using lightweight composite materials and a low-cost screen-printing process to inscribe the designed FSS pattern on an aramid sheet. Third, the proposed FSS pattern is an extension of the original idea in [START_REF] Barbagallo | Small periodicity FSS screens with enhanced bandwidth performance[END_REF] and propose to intertwine the resistive line strips on the four sides of the pattern to reinforce couplings and lowering the periodicity of the unit cell. This pattern restricts the meandering of strip lines for outmost ones, where it is the more relevant. A preliminary analysis of the candidate FSS structure was proposed in [START_REF] Merzaki | A compact absorbing FSS structure for antenna decoupling in the 5G 3.5GHz band[END_REF] and its operation and performance are further analyzed in details. The resistive surface of the decoupling device was fabricated using a screen-printing process with a home-made ink. This is the reason why recently proposed techniques based on graphene material [START_REF] Rinaldi | Graphenecoated honeycomb for broadband lightweight absorbers[END_REF]- [START_REF] Geng | A dynamically tunable microwave absorber based on graphene[END_REF] were not considered here. These solutions are not suitable to achieve the target sheet resistance value by the screen-printing process. The decoupling device is then obtained by stacking the screen-printed sheet on a composite substrate. The overall decoupling performance is optimized from numerical simulations and results are confirmed through experiments. The rest of this article is organized as follows. Section II deals with the design of the compact narrowband absorber made of an advanced FSS unit cell. It also describes its behavior when illuminated by a plane wave under oblique incidence. Section III details the fabrication process of the decoupling device. Section IV provides experimental evidence of its intrinsic absorbing performance. In Section V, the retained FSS structure is then shaped as a vertical wall for decoupling an array of stub-loaded antennas in the 3.4-3.6 GHz frequency band. Eventually, the experimental validation of the obtained decoupling device is presented. Section VI concludes the article.

II. FSS ABSORBER UNIT CELL DESIGN

An FSS absorber is composed of a periodic resistive pattern lying on a dielectric substrate with an electromagnetic reflector stacked on the opposite side. The contrast between the FSS absorber input impedance and the wave impedance determines the operating frequency band in which a minimum reflection coefficient is achieved. As described in [START_REF] Costa | Analysis and design of ultra thin electromagnetic absorbers comprising resistively loaded high impedance surfaces[END_REF], absorption takes place at the neighborhood of the resonance frequency of the stand-alone FSS pattern. The latter behaves as a simple resonant R, L, C circuit for a basic resistive pattern made of a patch, cross, or loop, where C and L stand for the mutual electric and magnetic coupling between adjacent cells, respectively, and R accounts for a resistance related to the pattern resistivity. The impedance at the FSS pattern interface is composed of the parallel association of the patterned surface and the impedance brought back from the grounded dielectric substrate [START_REF] Costa | Analysis and design of ultra thin electromagnetic absorbers comprising resistively loaded high impedance surfaces[END_REF]. Adjusting the thickness (which remains thin with regard to the wavelength), it is possible to obtain a narrowband absorber or a broadband absorber for a thicker dielectric. The maximum absorption occurs when the imaginary part of the two parallel impedances are opposite, the resistivity being adjusted to fulfill matching conditions.

The proposed FSS absorber was developed considering the original idea in [START_REF] Barbagallo | Small periodicity FSS screens with enhanced bandwidth performance[END_REF] that consisted in interdigitating adjacent cells. In this article, authors demonstrated that such pattern lowered the resonance frequency, increasing mutual coupling between contiguous cells. This suggests that intertwining neighboring cells should favors an increase of mutual coupling. While keeping the same approach, we therefore propose an evolution of the pattern, intertwining the resistive line strips on its four sides to reinforce couplings and lowering the periodicity of the unit cell. This pattern restricts the meandering of strip lines for outmost ones, where it is the more relevant. Furthermore, it was developed both to fulfill width resolution limitation (about 100 μm) and the range of attainable sheet resistance with home-made inks. Our design started from the elementary pattern in Fig. 1 resistive ribbons) as developed in [START_REF] Barbagallo | Small periodicity FSS screens with enhanced bandwidth performance[END_REF]. Our proposed TZZ inv (in series association of transformed Z-shaped and inverted Z-shaped resistive ribbons), structure initially presented in [START_REF] Merzaki | A compact absorbing FSS structure for antenna decoupling in the 5G 3.5GHz band[END_REF], is then presented in Fig. 1(b). It is first evaluated through a freestanding FSS simulation of an infinitely periodic and metallic (copper) FSS sheet in air in comparison with the initial ZZ inv pattern. Reducing the resonance frequency of an FSS structure leads to an intrinsic shift of the absorption frequency band toward lower frequencies when this resonant structure is transformed to FSS absorber, which is highlighted in the article. Note that both structures have the same periodicity (D = 5.9 mm). The geometry of the unit cell for the initial ZZ inv pattern and for the TZZ inv pattern is detailed in Fig. 1(c) and (d), respectively. Fig. 2 shows the physical view of the boundary conditions that has been ensured to simulate the FSS transmission coefficient of the free-standing pattern for both cases using the frequency solver of CST Microwave Studio. Simulations were carried out using two Floquet ports positioned on both sides of the z-axis. The boundary conditions are set at the unit cell level for x-, y-axis (shown in the green dashed rectangular parallelepiped). So, the structure is expanded periodically and infinitely in the xy plane. An open boundary condition is set along z-axis (blue squares) allowing incident wave to propagate through z-axis with minimum reflection before and after reaching the infinite structure. The FSS structure is illuminated by a normal incident plane wave source either in the transverse electric (TE) or in the transverse magnetic (TM) mode and located at the Floquet port 1. The transmission coefficient is therefore the transmitted power from the Floquet port 1 toward the Floquet port 2 located on the other side of the structure. Simulations of the transmission coefficient provide a better picture of the resonating nature of FSS patterns. The results of such simulations of both patterns are reported in Fig. 3.

The resonance frequency of the ZZ inv pattern, which is characterized by a null transmission coefficient, equals 10.44 GHz whereas the TZZ inv pattern exhibits a similar resonating behavior at a much lower frequency (6.17 GHz). Therefore, the decrease of the resonance frequency utterly comes from an increase of the mutual capacitance and inductance.

With the aim of developing a compact FSS absorber with an operating frequency centered at 3.5 GHz with a 200 MHzbandwidth, a parametric study of the dielectric permittivity is carried out to achieve the absorption peak at 3.5 GHz. Afterward, the sheet resistance optimization of the TZZ inv pattern is carried out to achieve the maximum absorption. The parametric study is performed with a 4 mm-thick dielectric substrate (λ/21.5 at 3.5 GHz).

Simulations are carried out by setting an unit cell boundary in the xy plane (at the four sides of the green dashed rectangular parallelepiped), and because of the presence of a reflector at the bottom of the absorber, a perfect electric boundary condition is set to feature a perfect ground plane (highlighted in red dash square, Fig. 4) in order to minimize the calculation time. Therefore, a single Floquet port is used as shown also in Fig. 4. The structure is illuminated by an incident plane wave source located at the Floquet port1 and the reflection coefficient is calculated at this port location.

Fig. 5 displays the TZZ inv simulations of 36 μm-thick copper ribbons (σ = 5.96 ×10 7 S/m) printed on a grounded 4 mmthick dielectric substrate. When the dielectric permittivity ε r increases, the absorption peak decreases smoothly toward low frequency. The absorption peak sets at 3.58 GHz with ε r = 3 and at 3.46 GHz with ε r = 4. Thus, to achieve the peak absorption at 3.5 GHz, the dielectric permittivity value could be expected to range between 3 and 4. To get the maximum absorption peak at 3.5 GHz, the dielectric permittivity value was set at 3.5 and new numerical simulations were performed by varying the sheet resistance value of the ribbons (see Fig. 6). With ribbons made of metallic conductor (copper with a sheet resistance R s = 0.0005 Ω/sq), the absorption peak reaches -3.98 dB at 3.47 GHz. When the sheet resistance value increases, the peak absorption level increases as the Ohmic loss in the TZZ inv pattern increases and the input impedance of the overall absorber matches with the free space impedance. This matching reaches its maximum at R s = 1.5 Ω/sq and the related absorption level achieves -18 dB at 3.47 GHz. Beyond 1.5 Ω/sq, the input impedance unmatched the free space impedance and the absorption peak level decreases.

Accepted manuscript

Based on the previous TZZ inv parametric studies, a compact absorbing FSS structure, called ATZZ inv was developed to operate in the 3.4-3.6 GHz frequency band. This absorber is characterized by a periodicity of only D = 5.9 mm (λ/14.5 at 3.5 GHz). The optimized sheet resistance of the ribbons is equal to 1.5 Ω/sq. The simulated substrate is a composite material made of glass fiber tissue and polyester resin. It is characterized by a dielectric permittivity ε r = 3.7, a loss tangent tanδ = 0.0094 and a thickness t = 4 mm. Fig. 7 provides an illustration of the simulated absorber whereas Fig. 8 shows the related simulated reflection coefficient at normal incidence over the 3-4 GHz frequency range. As a result, the absorption coefficient is better than -10 dB from 3.3 to 3.65 GHz. To assess the performance of the ATZZ inv FSS absorber, its simulated reflection coefficient was compared with that of the AZZ inv FSS absorber in the same configuration, i.e., using the same sheet resistance value of the resistive ribbons (R s = 1.5 Ω/sq) and the same dielectric substrate. At normal incidence, the AZZ inv cell based absorbing material reaches its maximum absorption at 6.5 GHz against 3.5 GHz for the ATTZ inv absorber. The gain in size is therefore of about 46%. Note that these simulations aim at comparing the position of the frequency absorption peak and not the absorption level. To compare the absorption levels of both structures, a sheet resistance optimization of the AZZ inv resistive ribbons needs to be performed.

Accepted manuscript

Fig. 9 illustrates the power loss density at 2 and 3.5 GHz of the ATZZ inv FSS absorber under the TE mode excitation. Power loss is mainly concentrated at upper and lower borders of contiguous FSS cells at 3.5 GHz compared with that at 2 GHz, which confirms the absorption effect at 3.5 GHz for the TE mode (electric field aligned with the y-axis). Under TM mode excitation, power loss is concentrated at the left and right borders of the unit cells (not shown here).

Such ATZZ inv absorber will be used in Section V to decouple antennas exposed to EM waves impinging at incidence angles different from normal incidence. Simulations are achieved using the boundary conditions of Fig. 4. By controlling a CST software predefined variable, the propagating wave direction is tuned to any required incident angle. Therefore, Fig. 10 illustrates the simulated reflection coefficient of the ATZZ inv FSS structure under four incident angles θ i in the TE mode (similar results, not shown here, are obtained in the TM mode). The absorption performance decreases under oblique incidence angle with an

III. FABRICATION PROCESS

The dielectric substrate is made of a fiber-reinforced matrix. The FSS pattern is printed on a dielectric sheet and stuck on the top face of the dielectric substrate.

On the one hand, the 4 mm-thick fabricated substrate (S2glass fiber tissue infused with polyester resin) was characterized at microwaves. Its dielectric permittivity ε r = 3.7 ± 0.5 and loss tangent tanδ = 9.4×10 -3 ± 0.5×10 -3 from 100 MHz to 1 GHz were retrieved from impedance measurements [START_REF] Sergolle | Structural composite laminate materials with low dielectric loss: Theoretical model towards dielectric characterization[END_REF]. The dielectric composite substrate is then associated with a 2 mmthick aluminum plate on its back-side (see Fig. 11).

On the other hand, the FSS pattern was printed using a screenprinting process based on a home-made ink (53.1% PF410 silver ink and 46.9% 6017SS carbon ink by weight) to get the required sheet resistance value of 1.5 Ω/sq (see Fig. 12). It is worth noting that the process resolution achieves an accuracy of 100 μm. The resistive FSS was patterned on a 100 µm-thick aramid sheet. The size of the aramid sheet, of the dielectric substrate and of the metallic plate equals 35 cm × 35 cm.

The absorbing performance of the fabricated absorber was measured using estimation of the radar cross section computed from S-parameters.

IV. MEASUREMENT METHOD AND EXPERIMENTAL RESULTS

Absorbing efficiency measurement is based on the retrieval of the radar cross-section (RCS) in an anechoic chamber (see Fig. 13). It is made possible through the retrieval of the ballistic backscattered wave from the target and by analyzing the scattering S-parameter as measured from TX and RX antennas in a quasi-monostatic configuration. The antennas were placed at a distance of 4 m from the sample. Identification of the target, considered as a point source scatter, was achieved using time gating of the inverse Fourier transform of the spectral response. The measurement of the absorbing FSS structure followed a twostep procedure. The first step consisted in measuring the RCS of the metallic plate. The second step consisted in measuring the FSS placed in front of such metallic plate. The contrast of RCS provides the absorbing performance of the FSS device, at normal incidence of the electromagnetic quasi-plane wave.

Fig. 14 shows the comparison of the simulated reflection coefficient of the two-dimensional infinite FSS absorber and the ratio between the RCS of the 35 cm × 35 cm FSS absorber and of the 35 cm × 35 cm metallic plate. This difference in the logarithmic scale is analogous to the magnitude of a reflection coefficient that can be directly compared with numerical simulation. The experimental result indicates that the realized absorbing FSS owns a maximum level of absorption at 3.43 GHz against 3.50 GHz for numerical simulation. The absorption bandwidth at -10 dB is 3.15-3.94 GHz, which remains wider than the simulated one (3.32-3.62 GHz). This difference is probably due to the difficulty in controlling the sheet resistance of the resistive ribbons all over the 35 cm × 35 cm sample size. Furthermore, the simulation is performed using periodic boundary conditions that do not take into account the finite size of the fabricated sample. As a reminder, TZZ inv absorption is sensitive to the ribbon sheet resistance value (see Fig. 6), a little variation in R s value modifying the behavior of the overall FSS structure.

Accepted manuscript

V. APPLICATION TO ANTENNA DECOUPLING

A. Stub Loaded Patch (SLP) Antenna Array

Utility of the FSS absorber for the reduction of the mutual coupling between antenna elements of an array is now assessed. An array of four patch antennas operating in the 5G 3.4-3.6 GHz frequency band is therefore considered. Each patch owns a set of four stubs attached symmetrically at the two opposite corners of their square patch. Hence, they are named SLP antennas. This four elements array is presented in Fig. 15. In the following, the four antennas are called Ant1, Ant2, Ant3, and Ant4 according to their labeling in Fig. 15(a). The physical dimensions of the antenna array are given in Table I. Its fabrication uses a composite dielectric substrate with a permittivity ε r = 2.4, a loss tangent of 0.0047 and a thickness t = 4 mm. Note that the spacing between antennas is only d 1 = 16 mm from edge to edge [see Fig. 15(a)].

Fig. 16 shows the magnitude of the reflection coefficient S 11 of Ant1 and the magnitude of scattering parameters from Ant1 to all other antennas, namely Ant2, Ant3, and Ant4. From the S 11 measurement, the operating frequency of the array lies between 3.38 and 3.62 GHz, S 11 being less than -10 dB over this frequency range. In that operating frequency band, the scattering coefficients to Ant2 (S 21 ) and Ant4 (S 41 ), i.e., in the direction of the E-plane cut of the radiation pattern, vary from -15 to -19 dB and from -25 to -40 dB, respectively. In the direction of the H-plane cut, the coupling coefficient S 31 varies from -18 to -22 dB. Specifications in terms of 5G MIMO antennas require that the coupling between each antenna remains at least lower than -20 dB. Therefore, to enhance our MIMO antenna system, the interposition of the proposed ATZZ inv FSS absorber is investigated to preclude the radiation coupling.

B. Decoupling Wall Absorber and Height Optimization

The decoupling device consists of the back-to-back association of two FSS absorbers of type ATZZ inv . They form a vertical wall between antennas to be decoupled. Fig. 17 shows the insertion of such decoupling wall in the E-plane of the

Accepted manuscript

TABLE II HEIGHT OF THE FSS DECOUPLING ABSORBER

antenna array where the coupling is significant. The thickness of the interposed device equals W = 8 mm, i.e., twice the thickness of a single absorber. Its length Le equals 150 mm. Its height was optimized to achieve a suitable tradeoff between significant decoupling achievements on the one hand and preservation of antenna radiation patterns on the other hand. Note that the separation from the decoupling device to each antenna is only 4 mm from edge to edge.

The optimization of the absorber height was achieved by varying the number of rows of the FSS design. The absorber performance was evaluated through the antenna input impedance, coupling rate, and radiation pattern impacts.

Fig. 18 shows the reflection coefficient of Ant1 for a two-rows (h 1 ), a three-rows (h 2 ) and a four-rows (h 3 ) configurations of unit cells. Table II gives the corresponding wall height values. Whatever the absorber configuration, the operating frequency of the array is slightly shifted to 3.4-3.63 GHz. Insertion of the wall absorber therefore enhances the matching performance, specifically above 3.45 GHz. More importantly, Fig. 19 highlights that the magnitude of the scattering parameter between Ant1 and Ant2 (S 21 ) is lower than -20 dB over the frequency range where Ant1 is matched for h 2 = 18 mm and h 3 = 24 mm. A very significant decoupling effect occurs beyond 3.5 GHz. At height h 3 , the decoupling wall operates as expected from 3.4 to 3.5 GHz where the decoupling is significantly improved compared with h 2 . Beyond the latter frequency, the coupling coefficient remains below -20 dB but the performance is much lower than that with height h 2 . We assume that diffraction by the edges of the absorbing wall causes this deterioration.

To ensure that decoupling reduction at height h 2 is lying on the absorption mechanism and not on diffraction effect, we performed additional numerical simulations removing either the FSS sheets or both sheets and dielectric substrates. Fig. 20 illustrates the transmission coefficients S 21 and S 31 where only the metallic plate remains, with addition of the substrates and finally with the complete back-to-back FSS decoupling wall. As far as S 21 is concerned, the FSS sheets allow a remarkable reduction of the transmission coefficient from 3.5 to 3.6 GHz compared with those of the metallic plate and dielectric substrate. This leads us to confirm that decoupling beyond 3.5 GHz is indeed due to the absorption effect. Moreover, the decoupling effect with regard to Ant3 is also demonstrated. The presence of the absorber wall shows a remarkable reduction of the coupling by an average of -4 dB from 3.4 to 3.5 GHz, once the FSS sheets are added onto the grounded dielectric substrates.

The height h 2 has been selected here as the optimized height of the FSS absorber, because of its performance achieved in the 3.5-3.6 GHz frequency band. Furthermore, its decoupling effect is linked to electromagnetic absorption and not to diffraction effect that may be detrimental to antenna radiation pattern. All transmission S-parameters from Ant1 to any other antennas are displayed in Fig. 21, with and without the h 2 -high decoupling wall. It appears that transmission coefficients remain below -20 dB for all parameters in the 3.43-3.6 GHz frequency range.

Fig. 22 highlights the current distribution generated by Ant1 with and without the h 2 -high FSS decoupling wall at 3.6 GHz. The current distribution of Ant1 is weakly modified by the presence of the absorber wall while current distributions induced in Ant 2 and Ant4 are significantly reduced.

C. Fabrication of the SLP Antennas Array

The four-elements antenna array was fabricated from a homemade 4 mm-thick composite laminate material. The dielectric substrate consists of Innegra fabric (woven polypropylene fibers, namely a 340 g/m 2 taffeta) infused with polyester resin. It exhibits a dielectric permittivity ε r = 2.4 ± 0.4 and a loss tangent tanδ = 4.7×10 -3 ± 0.5×10 -3 from 100 MHz to 1 GHz, values retrieved from impedance measurements. The conductive layers placed on both sides of the dielectric substrate are made of a 293 g/m 2 0.12 mm-thick copper woven fabric exhibiting a sheet resistance value R s = 6×10 -4 Ω/sq (measured from a standard 4-probe setup).

Accepted manuscript

The composite laminate material was fabricated by the vacuum infusion process. Five plies of Innegra fabric were stacked between two plies of copper fabric (single ground plane for the bottom; radiating elements and related feeding lines for the top). The whole was located on a glass slab. Placed under vacuum (-0.6 bar) with a plastic bag as shown in Fig. 23, the copper and Innegra fabrics were infused jointly with a liquid polyester resin (with 2% by weight of catalyst). After the resin polymerization at room temperature for 24 h followed by a postcure step at 40°C for 16 h, the composite laminate material, made of the dielectric substrate sandwiched between the two copper fabrics, was suitable for release. Afterward a three-axis machine was used to make the four microstrip antennas by subtractive method. The top copper fabric was therefore machined under 0.20 mm deep, removing out the unnecessary conductive layer. As shown in Fig. 24, the fabricated antenna array was characterized using a two-port VNA. Nonconnected antennas are left open.

Fig. 25 shows that the addition of the h 2 -high FSS decoupling wall does not alter the antenna array bandwidth (3.28-3.58 GHz) although a frequency shift is recorded with respect to numerical simulations.

On the one hand, the transmission coefficient from Ant1 to Ant2 (S 21 ) is mostly above -20 dB (-17 dB at maximum) over that frequency range if no decoupling device is present. On the other hand, the insertion of the h 2 -high FSS decoupling wall provides a significant decoupling of at least -21 dB and up to -25 dB over the same frequency range. The absence of a much higher decoupling performance as predicted at 3.6 GHz by numerical simulation (see Fig. 19) is supposedly due to the presence of the connectors at the onset of the antenna feeding lines. The decoupling performance toward Ant3 is enhanced of -4 dB from 3.5 to 3.55 GHz. Similarly, it is improved of about -4 dB at 3.4 GHz and -8 dB at 3.5 GHz toward Ant4.

In order to investigate the role of the connectors, E-plane radiation pattern of Ant1 was simulated and measured without the presence of the FSS decoupling wall at two frequencies: 3.5 and 3.6 GHz (see Fig. 26). At 3.5 GHz, the measured radiation pattern exhibits a minimum level at 30°of -10 dB, i.e., 6 dB below the simulated directivity in the same direction. This difference reaches -14 dB with regard to numerical simulation at 3.6 GHz. Furthermore, we have confirmed from another simulation (not shown here) that the presence of comparable metal pieces causes this radiation pattern alteration.

Accepted manuscript

Fig. 27 presents the E-plane and H-plane Ant1 radiation patterns at 3.5 GHz retrieved from simulations and measurements, with and without the h 2 -high FSS decoupling wall. The Accepted manuscript presence of the absorber causes a slight shift of the main lobe direction (4°) which is fully acceptable. Meanwhile, the H-plane radiation pattern has not been impacted by the presence of such decoupling wall. We may therefore conclude that this simple, compact and low-cost decoupling device enables to fulfill mutual decoupling requirements below -20 dB and, at the same time, fully preserves the antenna array performance.

VI. CONCLUSION

Mutual coupling between antennas of 5G base station is a critical issue with regard to beamforming capabilities. In this article, we have designed, fabricated, and investigated a specific composite absorber for mutual decoupling of an antenna array. This solution is based on a compact back-to-back FSS absorber wall of 8 mm-width and 18 mm-high above the planar antenna array.

The pattern of the unit-cell array was first optimized to reduce its size with regard to the wavelength, using intertwining cells. This led to the fabrication of an FSS absorber associating a composite laminate substrate fabricated by vacuum infusion process and a screen-printed sheet with an optimized sheet resistance value (R s = 1.5 Ω/sq). Intrinsic performance of the FSS absorber has been confirmed from experiments.

Then, the obtained absorbing structure was configured as a back-to-back FSS to form a decoupling wall for preventing mutual coupling of a planar SLP antenna array. The back-to-back FSS wall size was optimized through numerical simulations. Experimental validations are then performed to assess the mutual coupling reduction into the antenna array.

The achieved antenna array and decoupling device prototypes highlight the decoupling performance. Mutual coupling between antennas remains below -20 dB over the targeted frequency range (3.4-3.6 GHz), fulfilling this stringent requirement for so close antennas. At the same time, radiation pattern measurements also confirm that the proposed decoupling device do not introduced any significant alteration.

We believe that such a solution may be extended to various designs of antenna arrays in the frame of future deployment of 5G (and beyond) to reduce mutual coupling in base stations for MIMO antenna arrays.

Fig. 1 .

 1 Fig. 1. (a) Periodic basic FSS design and (b) proposed FSS design (D = 5.9 mm). (c) Dimensions of the basic unit cell. (d) FSS unit cell.

  (a) showing a freestanding FSS structure called ZZ inv (in series association of Z-shaped and inverted Z-shaped Accepted manuscript

Fig. 2 .

 2 Fig. 2. Physical representation of the TZZ inv unit cell free-standing simulation and the excitation modes of the plane wave.

Fig. 3 .

 3 Fig. 3. Simulated transmission coefficient T of the free-standing ZZ inv and TZZ inv patterns versus frequency.

Fig. 4 .

 4 Fig. 4. Physical representation of the TZZ inv absorber.

Fig. 5 .

 5 Fig. 5. Simulated reflection coefficient vs. dielectric permittivity values of a copper TZZ inv pattern printed on a 4 mm-thick dielectric substrate.

Fig. 6 .

 6 Fig. 6. Simulated reflection coefficient vs. sheet resistance values of conductive and resistive TZZ inv patterns printed on a 4 mm-thick dielectric substrate with (a) ε r = 3.5. Magnitude [dB] versus frequency and (b) Smith chart from 3 GHz to 4 GHz (b). The red star is located at 3.5 GHz.

Fig. 7 .

 7 Fig. 7. (a) ATZZ inv unit cell pattern and (b) ATZZ inv FSS absorber periodic pattern.

Fig. 8 .

 8 Fig. 8. Simulated reflection coefficient comparison vs. frequency of the AZZ inv and ATZZ inv absorbers at normal incidence.

Fig. 9 .

 9 Fig. 9. Power loss density of the ATZZ inv FSS absorber at (a) 3.5 GHz and (b) 2 GHz under TE mode excitation.

Fig. 10 .

 10 Fig. 10. Simulated reflection coefficient vs. frequency of the ATZZ inv FSS absorber under various incident angles θ i .

Fig. 12 .

 12 Fig. 12. Global view and detail of the ATZZ inv FSS absorber printed on a 100 μm-thick aramid sheet.

Fig. 13 .

 13 Fig. 13. RCS measurement setup.

Fig. 14 .

 14 Fig.14. Simulated reflection coefficient vs. frequency of the ATZZ inv FSS absorber compared with the difference between the RCS measurement of the metallic plate and that of the ATZZ inv FSS absorber placed in front of the metallic plate. Colored frames feature the simulated and measured absorption bandwidths (reflection coefficient lower than -10 dB).

Fig. 15 .

 15 Fig. 15. (a) Representation of the four-elements antenna array and (b) detailed view of a single antenna.

Fig. 16 .

 16 Fig. 16. Simulated S-parameters versus frequency of Ant1.

Fig. 17 .

 17 Fig. 17. Four-elements antenna array with the decoupling wall made of two ATZZ inv FSS absorbers.

Fig. 18 .

 18 Fig. 18. Simulated S 11 parameter versus frequency of Ant1 for various heights of decoupling wall (h 1 = 12 mm; h 2 = 18 mm; h 3 = 24 mm).

Fig. 19 .

 19 Fig. 19. Simulated S 21 parameter versus frequency of Ant1 for various heights of decoupling wall (h 1 = 12 mm; h 2 = 18 mm; h 3 = 24 mm).

Fig. 20 .

 20 Fig. 20. (a) Simulated S 21 and (b) S 31 parameters versus frequency with the metallic plate; with the metallic plate and the dielectric substrate; and with the complete FSS decoupling wall of height h 2 .

Fig. 21 .

 21 Fig. 21. Simulated S-parameters versus frequency with and without the h 2high FSS decoupling wall.

Fig. 22 .

 22 Fig. 22. Surface current distribution of the h 2 -high FSS decoupling wall at 3.6 GHz of Ant1 (a) with and (b) without absorber at 3.6 GHz.

Fig. 23 .

 23 Fig. 23. Vacuum infusion mounting to fabricate the home-made SLP antenna array.

Fig. 25 .

 25 Fig. 25. Measurement of the S-parameters from Ant1, with and without the h 2 -high FSS decoupling wall.

Fig. 26 .

 26 Fig. 26. Simulation and measurement of Ant1 radiation patterns in E-plane at 3.5 and 3.6 GHz without the FSS decoupling wall.

Fig. 27 .

 27 Fig. 27. Simulation and measurement of Ant1 radiation patterns in E-plane and H-plane at 3.5 GHz, with and without the h 2 -high FSS decoupling wall.

TABLE I PHYSICAL

 I DIMENSIONS OF THE SLP ANTENNA ARRAY

MERZAKI et al.: COMPACT DOUBLE-SIDED FSS ABSORBING WALL FOR DECOUPLING 5G ANTENNA ARRAYS

ACKNOWLEDGMENT

The authors would like to thanks J. Sol (IETR) for his technical assistance for microwave measurements and K. Desmars (SERIBASE Industrie) for his assistance in screen-printing process.

This work was supported in part by the French FUI STARCOM Project, in part the European Union through the European Regional Development Fund, in part by the Ministry of Higher Education and Research, in part by the Région Bretagne, in part by the Départements d'Ille-et-Vilaine and des Côtes d'Armor, and in part by Saint-Brieuc Armor Agglomération, through the CPER projects 2015-2020 MATECOM and SOPHIE.

Faissal Merzaki received the master's degree in electrical and electronic engineering in 2018 from Rennes1 University, Rennes, France, where he is currently working toward the Ph.D. degrees in electrical engineering.

His research interests include electromagnetic absorbers, specifically for the 5G and V\UHF application. At first, he traded international measuring equipment with Algeria Then after a stay abroad (England/Ireland). He was with the Séribase Industrie. During ten years, he developed the commercial activity in France and in Europe. Then, he joined the design office in 2015 in research, development and innovations.

Philippe

Thierry Levavasseur, photograph and biography not available at the time of publication.


Loreto Misuraca

Tatsuhito Matsuo

Aline Cisse

Josephine Loricco

Antonio Caliò

Jean-Marc Zanotti

Bruno Demé

Philippe Oger

Judith Peters

Jean- Marc Zanotti

High temperature molecular motions within a model protomembrane architecture

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

Introduction

The onset of first cellular life on Earth goes back to about 4 billion years 1 . Therefore, the origin and nature of the first cell membranes and their constituents are not fully determined. Because of the supposed lack of molecular complexity on the early planet, the first cell membranes were most likely composed of simple, single chain fatty acids 2 and fatty alcohols, raising questions as how they could withstand the very variable and extreme surrounding environment 3 of the oceanic hydrothermal vents or fields. One of the commonly accepted scenarios for the origin of life is indeed related to such places 4 . In particular, the membranes may have had to cope with high temperatures (typically around 85°C) and high hydrostatic pressures (up to a few hundred bars), and simultaneously accomplish elementary but crucial functions of cell membranes: divide the inside and outside into separate compartments, and assure a limited permeability combined with a certain stability and flexibility.

As proto-cells or proto-membranes are no longer accessible due to evolution, the only way to probe the properties of these hypothetical systems is to re-construct membranes from simple lipids with short chain lengths, presumably favoured by prebiotic synthesis [5][6][7] , and to expose them to the supposed environmental conditions. Here we concentrate on fatty acids and fatty alcohols and phospholipids, all with chain lengths of 10 carbons (named hereafter C10), which are among the shortest to self-assemble into lipid bilayers 8,9 .

Specifically, we study here several possibilities (see Figure 1): first, decanoic (or capric) acid, which is a C10 fatty acid known to form vesicles under fixed conditions 9 . Second, a possible architecture for protocell membranes, consisting of a 1:1 mixture of capric acid with a fatty alcohol of equal chain length (decanol). These short, single chain amphiphilic molecules were in fact the ones readily available in the prebiotic environment and have been found able to form stable vesicles 10 . In a recent work, we reported the remarkable properties of this membrane architecture which allow it to withstand extreme temperature conditions 11 . Indeed, such a combination of fatty acids and alcohols improves the vesicles' stability at room temperature and it induces a high-temperature conformational change at T ≥ 60 °C that leads to vesicle fusion with stable membranes at temperatures as high as 80 °C. In the following, this model membrane will be referred to as C10mix.

The third model is inspired by the work of Cario et al. 12 . These authors proposed an important role of long chain isoprenoid alkanes (referred as "apolar lipid") in the membrane mid-plane to further trigger its physico-chemical properties in relation to extreme conditions. In fact, Salvador-Castell et al. did find that they had a stabilizing role in model archaeal membranes 13 and Misuraca et al. reported an increased stability of lipid membranes in presence of alkanes 14 . Thus the C10mix was enriched with small percentages of eicosane, a linear alkane with 20 carbons (C20). The C10mix containing eicosane can be considered as the most complete prebiotic model that sums up the requested features of molecular abundance and need for a strategy against high temperature. To further compare the decanoic acid to a more evolved lipid, a phospholipid having two C10 acyl chains, but the much larger choline headgroup was chosen:

1,2-didecanoyl-sn-glycero-3phosphocholine (DCPC). The acyl chains are bound to the headgroup through a glycerol backbone, which can stiffen the structure. Not only the structural composition of lipidic membranes is significant for their behavior, but their dynamic properties have an equally important impact 15 . High temperature is known to increase the thermal energy and thus molecular motions within membranes, hence jeopardizing the stability of the system. Nonetheless, to be functional a cell membrane requires a specific state of fluidity and therefore high flexibility. Both together can be ensured by molecular dynamics within the membranes.

(a) (b) (c) (d)
Many kind of motions were identified within lipidic membranes, from small local motions up to collective motions of the whole system. A well adapted technique for their observation and analysis is incoherent neutron scattering, which probes essentially the motions of hydrogen atoms and of the molecular subgroups to which they are bound [16][17][18][19] . Neutron instruments give access to a broad range of energy resolutions, which are inversely propotional to the accessible time window, and thus allow different types of motions to be probed. In the present study, we combine elastic incoherent neutron scattering (EINS) experiments giving access to atomic mean square displacements (MSD) with quasi-elastic neutron scattering (QENS) experiments to characterize the kind of motion present in the sample and the associated amplitudes. To the best of our knowledge, such complete investigation on multilamellar vesicles of one type of lipids, of mixed lipids and alcohol or in presence of alkanes was never attempted before, mainly due to the lack of a suitable model to describe the dynamics. However, as molecular dynamics are so important for the understanding of the correct functioning in a certain state, we show here which parameters can be extracted, which ones are the most impacted by the structural differences and how we can interpret these variations.

Experimental section

Sample preparation

Decanoic acid, 1-decanol, perdeuterated eicosane, DCPC, bicine buffer and D2O were purchased from Sigma Aldrich (Merck). The samples were prepared by dissolving the decanoic acid and the decanol in a CHCl3 solution and DCPC in chloroform:methanol (2:1), followed by drying under a stream of nitrogen. The samples were then placed in a desiccator and left under vacuum overnight. The sample weights were checked at each step, to make sure that all the organic solvent was evaporated. For the samples consisting of decanoic acid : 1decanol (C10mix) the final ratio was 1:1.

The bicine buffer was prepared at a concentration of 0.2M by dissolving hydrogenated bicine in D2O, following previous protocols 10 . The buffer was filtered with a 0.2 μm millipore membrane before use. The dried organic solutions were then resuspended in buffer by vigorous vortexing for ≈ 1min, leading to the final milky-foamy solutions characteristic of the presence of large multi-lamellar vesicles (MLVs). The samples to be measured were then diluted in D2O buffer to a final concentration of 60 mg/ml (350 mM). The samples were placed in flat aluminium sample holders and sealed with indium. Their weights were checked before and after the experiment to control that no solution was lost.

Incoherent neutron scattering

The incoherent neutron scattering cross section of hydrogen exceeds by far that of all other atoms present in biological systems 20 and also that of its isotope deuterium. Therefore, contrast variation can be used to highlight specific parts of the samples, as in our case the lipids rather than the buffer or the alkanes, whose presence will only be accounted for by their impact on the membrane dynamics.

For our study, we used the spectrometer IN6-Sharp operated by the Laboratoire Léon Brillouin (Saclay, France) at the Institut Laue Langevin (ILL) in Grenoble, France. IN6-Sharp is a cold neutron time-of-flight 21 with an energy resolution of 70 µeV at 5.1 Å incident wavelength, giving access to a time window of about 10 ps. Elastic scans were performed from 278 to 355 K with a temperature step of 1 K and a measuring time of 2 minutes per point. QENS measurements were conducted on the same samples at the temperature points of 278, 293, 323 and 353 K for 2.5 hours each. Additionally, the bicine buffer was measured in the same conditions with the exception that the QENS scans were acquired for only 2 hours due to time limitations. An empty cell and the completely incoherent scatterer vanadium were measured at room temperature for correction and normalisation purposes. The empty cell and buffer data were subtracted from the sample data taking into account their different absorptions based on the correction formula of Paalman-Pings coefficients 22 . Vanadium data were used to normalise the elastic data and to determine the instrumental resolution for QENS data. Data reduction for EINS was carried out using the LAMP software 23 available at ILL and using the IGOR Pro software (WaveMetrics, Lake Oswego, OR, USA) for QENS.

Elastic incoherent neutron scattering

EINS measurements give access to the elastic incoherent structure function S(Q, 0 ± ΔE), where Q is the modulus of the difference between the incoming and outgoing neutron wavevectors and DE the energy resolution of the instrument.

The structure factor can be written, using the Gaussian approximation 24 , which assumes that the distribution of the atoms around their average position follows a Gaussian distribution,

𝑺(𝑸, 𝟎 ± ∆𝑬) ≈ 𝒆 ! 𝟏 𝟑 〈𝒖 𝟐 〉𝑸 𝟐 , (1)
Please do not adjust margins Please do not adjust margins where the <u 2 > are the MSD. As Q approaches zero, the approximation is strictly valid and it holds up to ⟨u 2 ⟩Q 2 ≈ 1 25 .

The average MSD can be obtained from the slope of the logarithm of the scattered intensities according to

〈𝒖 𝟐 〉 = -𝟑 𝝏(𝒍𝒏 𝑺(𝑸±∆𝑬)2 𝝏𝑸 𝟐 . ( 2
)
The Q-range used for fitting was on IN6-Sharp 0.45 to 1.4 Å -1 . It is possible to further determine the force constant or resilience <k> from the slope of the MSD 26 through

〈𝒌〉 = 𝟎.𝟎𝟎𝟏𝟑𝟖 𝒅〈𝒖 𝟐 〉 𝒅𝑻 . ( 3
)
Here <k> is expressed in Newton per meter when <u 2 > is given in Å squared and T is the absolute temperature.

Quasi-elastic neutron scattering

Whereas elastic scattering permits only averaged local motions to be seen, QENS is the tool to separate different kind of movements and to characterize them. For that, each spectrum can be fitted after the data reduction with a phenomenological equation containing terms which account for various diffusions:

𝑆(𝑄, 𝜔) = 𝐶(𝑄)+𝐴 ! (𝑄)𝛿(𝜔) + ∑ 𝐴 " (𝑄)𝐿 " (𝑄, 𝜔) " 1 × 𝑅(𝑄, 𝜔) + 𝐵(𝑄), (4)
G = 𝑫 𝑻 𝑸 𝟐 . ( 5
)
Sometimes, G is not crossing the origin for Q → 0, what is the signature of confinement for larger amplitudes. In this case, it is advised to fit the data only from a Qmin on, where the linear behaviour sets in.

A constant behaviour is typical for rotational diffusion characterized by the rotational diffusion coefficient DR:

G = 𝟐𝑫 𝑹 . (6)
Finally, jump-diffusion is typical for infinitely small, elementary jumps and a certain time spent by the protons in between, the residence time τ, in addition to the jump diffusion coefficient Djump. Here, G as function of Q 2 increases towards a constant limit, G∞ = 1/t:

G = 𝑫 𝒋𝒖𝒎𝒑 𝑸 𝟐 𝟏)𝑫 𝒋𝒖𝒎𝒑 𝑸 𝟐 𝝉 . ( 7
)
We fitted the experimental spectra in the Q-range from 0.37 to 2.02 Å -1 and within energy transfers w from -10 meV to 2 meV (see an example in Figure 2). Please do not adjust margins Please do not adjust margins

Results and discussion

Experiments were performed in August 2020 on IN6-Sharp (DOI: 10.5291/ILL-DATA.CRG-2728).

EINS

The integration of the elastic peak allowed the extraction of the MSD as a function of temperature through eq. ( 2) (see Figure 3). Generally speaking, lipids are very flexible molecules and our results are in close agreement with other examples in the literature (see reference [START_REF] Aoun | [END_REF] , but note that in this reference the use of a different factor for the MSD).

The slopes (<k>) of all four curves present very similar trends (see table 1), indicating that the four systems have similar stabilities. On the contrary, the MSD data highlights a few particularities for each one. As can be seen in Figure 3, the C10mix is the most flexible sample over the whole temperature range. However, we do not observe an impact of the conformational transition of the C10mix above 333 K on the local dynamics as seen by DLS 11 . By contrast, the capric acid seems to undergo a steep raise in MSD at around 305 K, indicative of a phase transition, which is somewhat higher than the one observed by DSC in previous studies (287K) 11,14 . We interpreted the phase at lower temperature as a collapsed phase, which melts around this temperature. The shift in temperature observed here can have various reasons such as different temperature ramps for the two techniques or differences in hydration [START_REF] Aoun | [END_REF] .

Please do not adjust margins

Please do not adjust margins characterized by a radius of confinement a, to be determined through a = p/Q0, where Q0 stands for the momentum transfer where jump diffusion sets in. In the present cases, this critical value of Q0 was the same for all samples and equal to 2.9 Å.

With respect to the global translational diffusion coefficient of the lipids within the membranes (see eq. ( 5)), as shown in Figure 4 as a function of temperature for the four samples, it seems to be very similar in almost all cases except for the sample C10mix with eicosane, which presents a much higher value at high temperature. Moreover, the G1 of this sample at 353 K shows a behaviour resembling more jump diffusion (increasing towards a limit and a constant value at the low-Q limit) than pure translational diffusion (linear behaviour), likely because of the presence of the alkanes which induces the confinement of the acyl chains, and limits the amplitudes of the motions. All values of the diffusion coefficients are in between 0.3 10 -5 cm 2 /s at 278K and 2 10 -5 cm 2 /s at 353K and lie below the value of bulk water, 2.3 10 -5 cm 2 /s at 298 K and 6.81 10 -5 cm 2 /s at 353K 30 , which represents somehow the upper limit for such diffusional processes.

The second type of motion corresponds to jump diffusion and allows the determination of the jump diffusion coefficient and residence time according to eq. ( 7). Figure 5 represents these quantities for all samples. However, this model does not describe very well the motions within C10mix + eicosane and DCPC, leading to huge error bars. We can speculate that 2, 3-site jump diffusion could be more appropriate here. The hierarchy of the residence times is consistent with the EINS results, e.g. C10mix < capric acid < C10mix + eicosane < DCPC. This quantity, representative for interactions between molecules, decreases with temperature as the lipids become more and more mobile and is highest for DCPC.

The last G3 permits to calculate the rotational correlation time τ3, defined as the inverse of G3, which is represented in figure S8 of the ESI. All values are very close and almost independent of temperature. These motions are the fastest with smallest amplitudes and thus almost not prone to the specific environment.

We further used a recently developed model to interpret the amplitudes A0 to A3 of eq. ( 4), the Matryoshka model 27 (see Figure 6). Although the Matryoshka model is dynamic, fitting the data with it, starting from values found in the literature, it allows to determine the impact of molecular dynamics on certain structural parameters. The amplitudes are globally fitted for each sample, but independently for each temperature point.

The Matryoshka model describes the two alkyl chains as one and the same cylinder. The modelled parameters are the following: R is the lateral diffusion radius of the whole lipid in the membrane plane. kforce represents a force constant or resilience characterising the relaxation of the lipid with respect to out-of-plane motions. The whole lipid can rotate, expressed through the distance head-lipid RH. Flip-flop motions of the head group are described by the angle a (we are using the definition of Pfeiffer et al. 16 and not the exchanges of lipids from one layer to the other, which correspond to much longer time scales) and the heads with radius bH can rotate around their own axis. R1 stands for the minimum radius of diffusion of the lipid tail, which according to the distance from the head increases with √𝒎, m being the index of the methylene group position on the tail and M the total number of these groups (here: M = 10 for all samples). H atoms inside the methylene Please do not adjust margins Please do not adjust margins and methyl groups can undergo jump-diffusion. The parameters needed to describe these motions are the distance H-C-H d between the two-sites and the probability of jump events f. The model depends further on experimental parameters known for lipids as the number z, which is the proportion of H atoms in the head (here z = 0.05 for capric acid, = 0.03 for C10mix with and without the eicosane, = 0.32 for DCPC) and on the fraction of H atoms, seen as immobile within the instrumental resolution, termed here pimm. For a validation of the model upon real data, see 31,32 .

In figure S9 of the ESI we show the fits of the amplitudes for all four temperatures of our samples. Globally, they work remarkably well, giving mainly rise to reduced c 2 values between 1 and 16. This is more satisfying as the model was first developed for phospholipids and successfully applied here for a first time to single chain fatty acids and alcohols, too.

The extracted fit parameters are presented in Figure 7. Some remarks are noteworthy: When z has a very small value (as here for capric acid and C10mix), it is almost impossible to model the head motions, which explains the huge error bars for parameters b and a for these samples. For DCPC, a was fixed to 30°, similarly to DMPC, as the head groups are the same 31 . In that case, the error bars are smaller and b is in the same order of magnitude as for DMPC 31 , which is consistent with identical head groups. The other quantity for which the fits cannot give reliable results is the tilt angle a. In fact, the head group of capric acid and C10mix is so small that its orientation does not make sense.

In contrast to DMPC 27 , R1 is close to 0 for almost every sample, meaning that the corresponding motion, the 2D diffusion of the effective tail, is negligible compared to the jump-diffusion of methyl groups. As it occurs for all systems with different head groups, this is probably related to the smaller chain length (M = 10, compared to M = 14 for DMPC), which would hinder the 2D diffusion of the tail. However, R1 seems slightly bigger for DCPC, which can be assigned to the lower phase transition temperature of DCPC so that these lipids are in a liquid crystalline state even at 293 K.

The jump-diffusion, represented by the f parameter, seems to be more likely to occur in DCPC, and much less in capric acid. This could be due to the larger head group and the liquid crystalline state of DCPC, allowing more motions within the chains, or to the fact that there are two tails in DCPC, compared to one in all other systems, and one effective tail considered in the Matryoshka model (as discussed by Bicout et al. 16 ). Eicosane does not seem to have a strong effect on this parameter, which is reasonable as it is supposed to be mainly in the midplane of the bilayers, and thus is not supposed to interact with the majority of methylene and methyl groups.

Overall, the fraction of mobile H atoms is similar for all samples, and is about 80 %. Notably, the hierarchy between the samples is close to the one encountered for the MSD in Figure 3, where the C10mix displays the highest mobile fraction, and DCPC the lowest. We also retrieve the lowering effect of eicosane, with a decrease in the mobile fraction, except for the highest temperature.

Concerning the motions of the whole lipid, we can first comment on the rotation, corresponding to the parameter RH.

While increasing with temperature, following the increase in thermal energy, the values of RH stand for all samples around 1.5 Å, which can be put in comparison with measurements by SANS and NMR. Indeed, they give an estimated value for the lipid height of 1.7 Å 11 . Here RH is close to the projection of the lipid height, which is then supposed to be slightly smaller. It is slightly smaller for DCPC than for the other samples, because the single chain amphiphiles will rotate along their longitudinal axis much more easily than DCPC 33 .

For the lateral diffusion radius and the in-out force constant, one recognises two groups, formed by DCPC and capric acid on one side and C10mix with or without eicosane on the other side.

We can conclude that in the latter case, in-plane motions of the lipids are more restricted, but the membrane is more flexible in the out-of-plane direction. The presence of eicosane can certainly explain such behaviour, but it seems already present in the C10mix alone. One might speculate that it favours the requirements of a protomembrane.

Please do not adjust margins

Please do not adjust margins

Conclusions

In summary, the exact composition of the lipid model membrane has a clear effect on molecular dynamics and structural parameters. EINS established the existence of a phase transition around 305 K for capric acid, likely corresponding to the melting of a collapsed phase already seen earlier 11 .

Otherwise we found the following hierarchy in the flexibilities with C10mix > capric acid > C10mix + eicosane > DCPC, although the stabilities (determined by the force constants, table 1) were very similar for all of them. It fosters the idea that the C10mix permits a high degree of adaptation to the environment due to the flexibility, which is however slightly decreased by the presence of the alkanes located in the midplane. As shown earlier 14 , the alkanes on the contrary increase the membrane's stability under extreme temperature and pressure conditions.

QENS measurements revealed that the C10mix sample containing eicosane presented a very high translational diffusion at the highest temperature (which is supposed to correspond to the native condition). In contrast to that, the jump diffusional coefficient, characteristic of smaller localised motions, of this sample evolved from a high value at low temperature to the same value as for capric acid and C10mix alone at high temperature. Thus such a small value seems mandatory for a good functioning. All these values are quite different from those of DCPC, which is a non-prebiotic phospholipid with a much bigger headgroup.

The amplitudes of the motions permitted further structural parameters to be obtained through a fitting procedure with some starting values taken from the literature, as observed through the dynamical changes. The two parameters b and a are not well defined for capric acid and C10mix due to their small head groups and are fixed for DCPC. Some parameters, such as RH, R1 and φ do not differ significantly among the samples. On the contrary, m, R‖ and kforce show respectively two groups of samples, assembling capric acid and DCPC on one side and C10mix and C10mix with eicosane on the other side. So one might speculate that these three parameters are the most important to guarantee the good functioning of a membrane built from short-chain lipids and alcohol, with or without alkanes in the midplane. We retain a higher mobile fraction for these samples, especially at higher temperature, paired with a smaller lateral diffusion radius and a lower in-out force constant, indicating the motional direction perpendicular to the surface being favoured over lateral motions.

Figure 1 :

 1 Figure 1: Four models for protomembranes: (a) Capric acid, (b) Capric acid + decanol (C10mix), (c) C10mix + eicosane, (d) DCPC.

  where w is the energy transfer in multiples of ħ, C(Q) is a scaling factor to adjust the model and the absolute data, and A0(Q) δ(w) is the Q-dependent elastic incoherent structure factor (EISF). The functions Li(Q,w) are Lorentzian functions of half-width at half-maximum (HWHM) GI, each describing a certain type of diffusive motion in the lipids. Their amplitudes Ai(Q) are also Qdependent. R(Q,w) refers to the instrumental resolution which is obtained from the vanadium measurements. B(Q) is a flat background. The sum can theoretically be extended over an infinity of contributions, but real data allow only to resolve a few of them. Careful checks determined that 3 Lorentzians were sufficient to reasonably fit the data here without overfitting them. They represent slow (e.g. two-dimensional diffusion within the membrane plane and one-dimensional diffusion of the whole lipid perpendicular to the plane), intermediate (e.g. rotational diffusion of the whole lipid, flip-flop and rotational motions of the headgroup) and fast motions (e.g. twodimensional diffusion of the tail group and jump-diffusion of hydrogen atoms in methyl groups), as classified in the Matryoshka model recently proposed for a very complete description of motions within a membrane27 .The behaviour of G as function of Q 2 is indicative for the type of movement 28 : a linear behaviour represents translational diffusion characterized by the diffusion coefficient DT:

Figure 2 :

 2 Figure 2: Examples of fit results for the spectra of (a) capric acid and (b) C10mix, both at Q = 1.38 Å -1 and T = 293 K. The grey circles are the data points with errors, the black line is the total fit line. The green line represents the elastic contribution convoluted by the resolution function. The magenta, orange and blue curves are the Lorentzian functions corresponding to slow, intermediate and fast motions, respectively, convoluted by the resolution.

Figure 3 :

 3 Figure 3: MSD for MLVs of the four samples as a function of temperature.

Figure 4 :

 4 Figure 4: Global translational diffusion coefficient of all four samples as a function of temperature.

Figure 6 :

 6 Figure 6: Representation of a lipid molecule within the Matryoshka model and fit parameters.

Figure 5 :

 5 Figure 5: The jump diffusion coefficient and the residence time for all samples as a function of temperature extracted from G2.

Table 1 :

 1 The force constants of the four samples.In presence of eicosane, the C10mix becomes less flexible, as expected due to the hindrance of motions by an obstacle in the midplane of the bilayer, but as mentioned above, it maintains a stability similar to the C10mix as indicated by the slope values. Finally, DCPC presents the lowest MSD, only coinciding with those of capric acid at the lowest temperatures. The bigger headgroups, and the connection of the two C10 chains to the glycerol backbone could indeed hinder motions at the level of the heads or within the hydrophobic region.

		<k> (N/m)	D<k> (N/m)
	C10mix	0.102	0.003
	C10mix + eicosane	0.092	0.002
	DCPC	0.096	0.002
	Capric acid (< 305 K)	0.085	0.005
	Capric acid (>305 K)	0.124	0.007

K, in accordance with the phase transition found from EINS. The C10mix seems again the most flexible, but very similar to its counterpart with eicosane. The amplitudes A0 to A3 will be treated in more detail with the Matryoshka model 27 to extract also structural information (see below).

Second, the HWHM GI are representative of three identified types of motions (see Figures S5 to S7 in the ESI). It can be deduced by the different energy ranges they cover:

This journal is © The Royal Society of Chemistry 20xx J. Name., 2013, 00, 1-3 | 3

This journal is © The Royal Society of Chemistry 20xx J. Name., 2013, 00, 1-3 | 7

Acknowledgements

This work was funded by the French National Research Agency program ANR 17-CE11-0012-01 to P.O. and J.P. L.M. was supported by a scholarship from the Institut Laue-Langevin (ILL) PhD program. The authors thank ILL for neutron beamtime on IN6-Sharp (DOI: 10.5291/ILL-DATA.CRG-2728).

Figure 7: Fit parameters obtained for all samples as a function of temperature. Error bars are within symbols if not shown.

Conflicts of interest

There are no conflicts to declare.


Johann Priou
email: johann.priou@inria.fr

Szymon Greś

Matthieu Perrault

Laurent Guerineau

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

Keywords: Operational Modal Analysis, uncertainty quantification, stabilization diagram, automated interpretation, hierarchical clustering, subspace methods

The interpretation of stabilization diagrams is a classical task in operational modal analysis, and has the goal to obtain the set of physical modal parameters from estimates at the different model orders of the diagram. The diagrams are contaminated by spurious modes that appear due to the unknown (non-white) ambient excitation and sensor noise, as well as possible over-modelling. Under the premise that spurious modes will vary and physical modes will remain quite constant at different model orders, the focus is to retrieve the physical modes that constitute the identified model, while rejecting the non-physical, spurious modes. Over the last decade, extensive research has been devoted for developing automated strategies facilitating their interpretation. To this end, the interpretation is in principle disconnected from the identification method and boils down to three stages i.e., clearing the diagram from the spurious mode estimates, aggregating the modal parameter estimates in modal alignments and the final parameter choice. Besides the point estimates of the modal parameters, also their confidence bounds are available with some identification methods, such as subspace identification. These uncertainties provide useful information for an automated interpretation of the stabilization diagrams. First, modes with high uncertainty are most likely non-physical modes. Second, the confidence bounds provide a natural threshold for the automated extraction of modal alignments, avoiding the requirement of a deterministic threshold regarding the allowable variation within an alignment. In this paper, a strategy is presented for the automated mode extraction considering their uncertainties, based on clustering a statistical distance measures between the modes. The relevance of the uncertainty consideration in the automated extraction will be demonstrated on vibration data from two bridges.

INTRODUCTION

Operational modal analysis is a fundamental task in engineering practice. Assuming vibrating structures with a linear time-invariant behaviour, linear system identification techniques such as subspace-based system identification [START_REF] Peeters | Reference-based stochastic subspace identification for outputonly modal analysis[END_REF][START_REF] Döhler | Fast multi-order computation of system matrices in subspacebased system identification[END_REF] can be used, where the system matrices of a state-space model are estimated before retrieving the modal parameters from their eigenstructure. The model order of the system (corresponding to the number of modes in the data) is in general unknown, and moreover non-physical modes due to colored noise are usually present besides the physical structural modes of interest. Under the premise that spurious modes will vary and physical modes will remain quite constant at different model orders, the identification is repeated at different model orders using efficient algorithms [START_REF] Döhler | Fast multi-order computation of system matrices in subspacebased system identification[END_REF] to obtain the stabilization diagram. Then the focus is to retrieve the physical modes that constitute the identified model, while rejecting the non-physical, spurious modes. With subspace-based methods, not only the modal parameter estimates but also their uncertainty bounds can be obtained [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][START_REF] Mellinger | Variance estimation of modal parameters from output-only and input/output subspace-based system identification[END_REF][START_REF] Greś | Uncertainty quantification of modal parameter estimates obtained from subspace identification: An experimental validation on a laboratory test of a large-scale wind turbine blade[END_REF], which are valuable information for clearing the diagrams of modes with high uncertainty and their further processing.

Over the last decade, extensive research has been devoted for developing automated strategies facilitating the interpretation of the stabilization diagrams, e.g. [START_REF] Magalhães | Online automatic identification of the modal parameters of a long span arch bridge[END_REF][START_REF] Reynders | Fully automated (operational) modal analysis[END_REF][START_REF] Yaghoubi | Automated modal parameter estimation using correlation analysis and bootstrap sampling[END_REF], where an overview is given in [START_REF] Zini | A quality-based automated procedure for operational modal analysis[END_REF]. Such strategies usually comprise three stages i.e., clearing the diagram from the spurious mode estimates, aggregating the modal parameter estimates in modal alignments and the final parameter choice. A popular tool for both the clearing and the aggregating stages is clustering based on a distance measure between the modes of the stabilization diagram. A drawback of clustering in this context is that either the number of clusters i.e., the number of modes, is usually assumed to be known a priori, or that thresholds on the allowable distance of the elements within a cluster need to be defined.

The focus of this work is to alleviate these issues to achieve less dependence on user-defined thresholds. First, a distance measure for clustering is proposed that is based on the modal parameters and their uncertainties. In this way, the natural variation of the modes between different model orders due to the intrinsic uncertainties is considered. Thresholds for the distances can then be defined based on a chosen confidence level on statistical grounds, instead of defining them in absolute terms on allowable changes in the modal parameters that could be quite different for different noise levels of the data. Second, to further reduce dependence on the user-defined thresholds, an automatic post-processing of the clusters with a merging and a separation step is introduced based on the cluster distances, where too tight or too loose thresholds are alleviated for the final clusters. This paper is organized as follows. In Section 2., the basics of the subspace-based system identification and uncertainty quantification are recalled. In Section 3., the clustering strategy is developed for the automated interpretation of the stabilization diagrams and applied to vibration data of two bridges in Section 4., before concluding the paper in Section 5..

SUBSPACE-BASED SYSTEM IDENTIFICATION AND UNCERTAINTY QUANTIFICATION

System identification

Assume that the vibration behavior of the investigated structure can be modelled by a linear timeinvariant system, then its dynamics can be described by the discrete-time state space model

x k+1 = Ax k + w k y k = Cx k + v k , ( 1
)
where A is the state transition matrix, C is the output matrix, and k is the integer time step corresponding to the system at time t = k∆t, where ∆t is the sampling rate. Vector y k ∈ R r contains the measured outputs (such as accelerations, velocities, displacements, strains), and x k ∈ R n is the state vector. The state noise w k ∈ R n is related to the unknown ambient excitation, and vector v k ∈ R r is the output noise.

The modal parameters are related to the eigenvalues and eigenvectors (λ i , φ i ), i = 1, . . . , n, of A and to C by

µ i = log(λ i ) ∆t , f i = |µ i | 2π , ζ i = -Re(µ i ) |µ i | , ϕ i = Cφ i , (2)
where µ i is an eigenvalue of the corresponding continuous-time system, f i is the natural frequency, ζ i is the damping ratio and ϕ i is the mode shape at the output coordinates.

To estimate the system matrices A and C from the output data y k of length N , k = 1, . . . , N , and consequently the modal parameters in (2), the reference-based covariance-driven subspace algorithm [START_REF] Peeters | Reference-based stochastic subspace identification for outputonly modal analysis[END_REF][START_REF] Döhler | Fast multi-order computation of system matrices in subspacebased system identification[END_REF] is used. The output covariance estimates with respect to a subset of reference sensors or projection channels are computed as Ri = 1

N N k=1 y k+i y (ref)T k
. Arranged in block Hankel format, their theoretical values satisfy the decomposition

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      = O p+1 C q , where O p+1 =      C CA . . . CA p      (3)
with the observability and stochastic controllability matrices O p+1 and C q , and where p, and q are time lag parameters. With the singular value decomposition (SVD)

Ĥp+1,q = U 1 U 0 S 1 0 0 S 0 V T 1 V T 0 , an estimate Ôp+1 = U 1 S 1/2 1
of the observability matrix is obtained for the model order corresponding to the truncation of the SVD. Then, the output matrix C is estimated from the first block row of the observability matrix, and the state transition matrix is estimated from the shift-invariance property  = ( Ô↑ p+1 ) † Ô↓ p+1 where Ô↑ p+1 and Ô↓ p+1 are the observability matrix estimate without the last and first block row, respectively. Ultimately, the modal parameters are obtained from the eigenvalues and eigenvectors of  and from Ĉ as in (2).

Uncertainty quantification

The computation of the modal parameter covariance results from the propagation of the sample covariance on the Hankel matrix estimate Ĥp+1,q through all steps of the modal identification algorithm. This sample covariance reflect in particular the unknown inputs due to non-measurable excitation sources and the sensor noise, and they contribute in a non-trivial way to the covariance of the modal parameter estimates. The estimation of the sample covariance of the Hankel matrix is straightforward and is done by separating the dataset for k = 1, . . . , N into blocks. The propagation to the modal parameter estimates is then based on the delta method [START_REF] Casella | Statistical inference[END_REF], where the analytical sensitivity matrices are obtained from perturbation theory [START_REF] Reynders | Uncertainty bounds on modal parameters obtained from stochastic subspace identification[END_REF][START_REF] Mellinger | Variance estimation of modal parameters from output-only and input/output subspace-based system identification[END_REF].

Let ∆X be a first-order perturbation of a matrix-valued variable X, like X = H p+1,q . Then, for a function Y = f (X) it holds vec(∆Y ) ≈ J Y,X vec(∆X), where J Y,X is the derivative that can be obtained by analytically perturbing the functional relationship between X and Y , and where vec(•) denotes the column stacking vectorization operator. Subsequently, covariance expressions for the estimates satisfy

cov(vec( Ŷ )) ≈ J Y,X cov(vec( X))J T Y,X . (4)
In this context, the modal frequencies and damping ratios satisfy

∆f i ≈ J f i ,H vec(∆H), ∆ζ i ≈ J ζ i ,H vec(∆H), (5)
where the analytical sensitivities are derived in detail in [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 modal parameter covariance follows from (4) based on the Hankel matrix sample covariance and estimates of the respective sensitivities.

AUTOMATED INTERPRETATION OF STABILIZATION DIAGRAMS

Clustering techniques are well-suited for the automated interpretation of stabilization diagrams [START_REF] Reynders | Fully automated (operational) modal analysis[END_REF][START_REF] Pereira | Modal tracking under large environmental influence[END_REF], where data (the modes in the stabilization diagram) are regrouped according to similar properties into clusters (the modal alignments) with respect to a chosen distance measure. If the distance between two modes is "small", they are grouped together in one cluster. The clusters should contain modes with stable properties over several model orders. This section will first define appropriate distance measures, and in a second part will describe the clustering process to obtain alignments.

Definition of distance

Amongst several distances for clustering [START_REF] Reynders | Fully automated (operational) modal analysis[END_REF], a basic distance is based on the frequency difference and the MAC between two modes i and j in the stabilization diagram, defined as [7]

D c ij = |f i -f j | max(f i , f j ) + (1 -M AC ij ). (6)
Here the absence of the damping ratio can be noticed, which is often not used for clustering due to its high variability. The distance is normed in the sense that both the (normed) frequency and the MAC are in the interval [0, 1], leading to

D ij ∈ [0, 2]
where D ij indicates a perfect match between the modes.

Despite its high variability, the damping ratio is an important information to be considered in the distance measure, which can be achieved by appropriate weighting in

D n ij = α 1 |f i -f j | max(f i , f j ) + α 2 |ζ i -ζ j | max(ζ i , ζ j ) + α 3 (1 -M AC ij ), (7)
with the weighting factors α l ∈ [0, 1] with α l = 1 the distance can be interpreted as a change in percent between two modes. In the classical distance (6) a 50% weighting is chosen by default. Here a weighting is chosen as α = [50%, 5%, 45%] between the modal parameters, reflecting the high uncertainty of the damping ratio. The threshold of the distance for elements in the same cluster is chosen ad-hoc as D lim = 0.9, corresponding to 10% allowable variation in each of the weighted modal parameter components.

The normed distance [START_REF] Magalhães | Online automatic identification of the modal parameters of a long span arch bridge[END_REF] is based on the (deterministic) relative difference between the modal parameters. Since the modal parameter uncertainties are available from their estimation as outlined in Section 2.2., this uncertainty can be used as a natural weighting in the computation of the distance. Considering the uncertainties of the frequencies and damping ratios, the first part of the distance can be evaluated statistically by

S ij = f i -f j d i -d j T Σ -1 ij f i -f j d i -d j (8)
where

Σ ij = E ∆(f i -f j ) ∆(d i -d j ) ∆(f i -f j ) ∆(d i -d j ) T ,
which can be evaluated based on Equations ( 4) and [START_REF] Mellinger | Variance estimation of modal parameters from output-only and input/output subspace-based system identification[END_REF]. Since the modal parameter estimates are asymptotically Gaussian distributed [START_REF] Mellinger | Variance estimation of modal parameters from output-only and input/output subspace-based system identification[END_REF], the variable S ij is asymptotically χ 2 distributed with two degrees of freedom. Its non-centrality parameter is different from zero if the theoretical values of the frequencies and damping ratios are different. Hence, a known value for a target percentile of the distribution is available for verifying the assumption that the theoretical values of the frequency and damping ratios are equal, thus belonging to the same cluster. When targeting a 99% percentile, the corresponding threshold obtained from the χ 2 distribution is around 9.2, so modes with S ij < 9.2 can be considered to be a member of the same cluster. A criterion for the M AC can also be added on statistical grounds based on [START_REF] Greś | Uncertainty quantification of the Modal Assurance Criterion in Operational Modal Analysis[END_REF]; however, for ease of presentation its deterministic format is kept in the following, defining the combined distance, also called statistical distance,

D s ij = S ij + α(1 -M AC ij ). (9)
with an appropriate weighting α. This weighting factor can be chosen as follows. Considering a weighting of 45% of the MAC to the distance measure [START_REF] Yaghoubi | Automated modal parameter estimation using correlation analysis and bootstrap sampling[END_REF] as in [START_REF] Magalhães | Online automatic identification of the modal parameters of a long span arch bridge[END_REF], and the remaining 55% to S ij having a threshold of 9.2, a total threshold for D s ij in ( 9) is set to D s lim = 16. Then, assuming a criterion for the allowable change of the MAC in the order of M AC ij > 0.9, this corresponds to α = 70.

Clustering

There are many methods for clustering in the literature; amongst the most known are k-means and hierarchical clustering [START_REF] Duda | Pattern Classification (2nd)[END_REF] which are both in use for automated modal analysis. The k-means clustering requires the number of clusters as an input, which is hard to know beforehand in modal analysis. However, it is mostly used to distinguish spurious modes from physical ones [START_REF] Reynders | Fully automated (operational) modal analysis[END_REF][START_REF] Neu | Fully Automated Operational Modal Analysis using multi-stage clustering[END_REF]. Hierarchical linkage clustering is based on a tree format where elements are clustered according to their distance, starting with the smallest one and stopping when a maximal distance D lim is reached. In the first step, each element defines a cluster. In an iterative process, the two clusters with the shortest distance are combined into a larger cluster, where the distance between two clusters is defined in different ways. Setting a limit on the distances between the modes instead of defining the number of modes seems to be more natural for clustering the stabilization diagram; hence the hierarchical linkage clustering is chosen in the following.

Linkage

The linkage is the choice of how the elements (i.e., the modes) form a cluster. Two popular ways of linkage are complete and single. Single linkage regroups elements together by the minimum distance, i.e., the distance between two clusters is defined by the minimum of the distances between their respective elements when combining them. It can also be called nearest neighbor. On the other side, the complete linkage regroups the clusters by the maximum distance of their elements, also known as farthest neighbor. The direct consequence is that the complete linkage is more strict as the maximum distance of elements should be under a specified value, while the single linkage allows more flexibility. For the first clustering step, complete linkage is chosen to ensure that the distances between all cluster elements are below the respective threshold D lim as defined previously. Resulting clusters with too few elements are deleted. In a second step, a post-processing of the retained clusters is carried out. To avoid modal alignments that are split over the model orders due to too low thresholds, the clusters are checked with the single linkage distance for merging. Furthermore, to avoid the contrary effect of over-agglomerated clusters with several modes at single model orders due to too high thresholds, a separation step is introduced.

Post-processing

As presented before, the clustering depends on the threshold D lim which is an a priori choice to a greater or lesser extent -greater when considering deterministic distance measures and lesser when considering statistical distance measures. If the threshold is too low physical modes will be missing; if it is too high different modes will be agglomerated in the same cluster. For less impact of the choice of this parameter several correction steps have been developed.

Merge

The first correction step is merging. The complete linkage for the initial clustering is strict and may tends to divide single mode alignments into different clusters. To correct this, a near neighbor merging was implemented comparing close neighbors between them based on the single linkage approach. The minimum distance between two close clusters is computed, then if this distance is below the limit D lim the two cluster are merged.

Separation

The second step is a separation. A global mode cannot have different estimates at the same model order, yet this can happen when the threshold is too big, or when noise modes with similar properties are nearby. The purpose is then to clean some outlier points that could have been mixed with a real mode, and to possibly separate two modes that are close. To clean these clusters, a reference "mono-order" subset is defined by the elements of the cluster with unique model orders, i.e. at the corresponding model orders there is only one mode estimate in the cluster. For the remaining model orders with multiple estimates the closest element to the reference sub-set is kept, and the others deleted if too close.

The final step is a size verification of the found cluster. Every cluster whose number of elements is under a stability criterion nb stab is deleted. The criterion nb stab is chosen as 20% of the maximum order computed. The modal parameters are then determined from the means of the cluster elements when using the classical normed distance, or from the weighted mean of the elements, where the weights are related to the modal parameter uncertainties at the different model orders [START_REF] Döhler | Variance computation of modal parameter estimates from UPC subspace identification[END_REF]. while There are model orders with several modes do 7:

Algorithm
Compare the different modes at each order to the mode set at orders with single modes 8:

Keep the closest and delete the other 9: end while 10: end for

APPLICATIONS

Vibration data of two bridges, namely S101 Bridge [START_REF] Döhler | Structural health monitoring with statistical methods during progressive damage test of S101 Bridge[END_REF] and Z24 Bridge [START_REF] Reynders | Continuous vibration monitoring and progressive damage testing on the Z24 bridge[END_REF][START_REF] Döhler | Uncertainty quantification for modal parameters from stochastic subspace identification on multi-setup measurements[END_REF] are used to apply the developed cluster algorithm 1 using either of the two distances defined in [START_REF] Magalhães | Online automatic identification of the modal parameters of a long span arch bridge[END_REF] and in [START_REF] Yaghoubi | Automated modal parameter estimation using correlation analysis and bootstrap sampling[END_REF], denoted by normed distance D n and statistical distance D s , respectively.

S101 Bridge

One dataset of the reference state of S101 Bridge before damage [START_REF] Döhler | Structural health monitoring with statistical methods during progressive damage test of S101 Bridge[END_REF] has been chosen for analysis, containing measurements from 45 sensors during 5.5 min. The respective stabilization diagram with the uncertainty bounds on the frequencies is shown in Figure 1, where the modal alignments are already quite clear.

In Figure 2, the clustering results with the developed algorithm are shown based on the classical normed distance (top) and the statistical distance (bottom). The five modes that were previously identified in [START_REF] Döhler | Structural health monitoring with statistical methods during progressive damage test of S101 Bridge[END_REF] are shown as vertical lines in Figure 2. They are clearly present with both distances. This validates the use of statistical distance.

The uncertainty bounds of the modal parameters can also be used as hard criterion to delete very uncertain modes.With a threshold on the coefficient of variation of the frequencies (standard deviation of the frequency divided by the frequency), the stabilization diagram becomes clearer as shown in Figure 3 in Figure 4, where they do not seem to be impacted by the application of the uncertainty threshold, while allowing as faster computation of the clustering.

The clustering algorithm is efficient to automatically retrieve stable alignments; however some alignments of spurious modes are also stable but tend to vary between different datasets. In order to retrieve the modes that are stable over different datasets, the same clustering technique can be applied on the set of modes that is obtained from each of the datasets in the spirit of "cross-validation" of the modes. Assuming only little change of the physical modes between the analyzed datasets, the same threshold for their distances as in the analysis of a single stabilization diagram can be used.

The results of this second clustering on the modal parameters obtained from 21 successive datasets in the reference state of S101 Bridge is shown in Figure 5 well recognized and other modes with higher variability are dismissed. To account for higher mode variability between different datasets (as e.g. in tracking applications [START_REF] Pereira | Modal tracking under large environmental influence[END_REF]), the threshold is usually increased. Results with the doubled threshold are shown in Figure 5 (right), where more stable modes are found. An in-depth analysis of these additional modes showed that the modes at 1 Hz and 18 Hz are most probably non-physical, while the other modes are indeed physical modes with lower excitation.

Z24 Bridge

The second example is the Z24 Bridge [START_REF] Reynders | Continuous vibration monitoring and progressive damage testing on the Z24 bridge[END_REF][START_REF] Döhler | Uncertainty quantification for modal parameters from stochastic subspace identification on multi-setup measurements[END_REF], where one dataset of the multi-setup measurements has been used as in [START_REF] Döhler | Efficient multi-order uncertainty computation for stochastic subspace identification[END_REF]. Measurements are more noisy, resulting in higher uncertainties of the modal parameters and more spurious modes, as shown in Figure 6, making the automated analysis more difficult. Nevertheless, the same thresholds as in the previous example are used for the clustering, and the impact of the developed steps for post-processing (merge and separate, see Algorithm 1) as well as the impact of the statistical distance measure are analyzed. In Figure 7 the result of the clustering is shown without any of the post-processing steps described in Section 3.2.2.. The black vertical lines represent the reference modes identified in [START_REF] Döhler | Uncertainty quantification for modal parameters from stochastic subspace identification on multi-setup measurements[END_REF]. Indeed, some of the alignments are split due to the more noisy data that would possibly require more relaxed thresholds. After the post-processing steps, the final clusters correspond well to the structural modes in Figure 8. However, there is one exception, namely the mode at 17 Hz that is very lowly excited and shows high variability. It is missing when using the classical normed distance; however the statistical distance considers its higher uncertainty intrinsically where it is well recognized. This shows that the statistical distance is more robust in a noisy situation.

It should be noted that it the effect of cleaning the diagram was tested by putting a threshold on the coefficient of variation of the frequencies, similarly as in Figures 2 and3. First, since data is noisier, the threshold of 1.5% needed to be increased in order to keep modes with higher estimation uncertainty. Second, the results did not show improvement of the clustering, as already seen in Figure 3; moreover the mode at 17 Hz (which has a higher variation) tends to get lost if the threshold is not set very carefully. This shows that clustering strategy does not require a previous cleaning of the diagram, whereas such a cleaning may be helpful for a manual visual inspection of the diagram.

CONCLUSIONS

The proposed statistical distance for clustering of stabilization diagrams has been validated on the S101 Bridge data. Then an application on the Z24 Bridge showed that in difficult noise condition the statistical approach allows to find modal alignments even of very noisy modes due to its intrinsic consideration of the uncertainties, leading to a more robust algorithm. The examples also showed that clearing the diagrams with a hard criterion on modal uncertainty bounds may not be necessary, since the clustering algorithm takes them into account directly, and may even be counterproductive. Furthermore, it has been shown on the S101 Bridge data that the clustering approach is a means for cross-validation of modes from several datasets.

Figure 1 :

 1 Figure 1: Full stabilization diagram of the S101 Bridge

Figure 2 :

 2 Figure 2: Hierarchical clustering from full stabilization diagram of the S101 Bridge (Figure 1) based on normed distance D n (top) and statistical distance D s (bottom)

Figure 3 :

 3 Figure 3: Cleaned stabilization diagram of the S101 Bridge with threshold on the coefficient of variation of frequencies of 1.5%

Figure 4 :

 4 Figure 4: Hierarchical clustering from cleaned stabilization diagram of S101 Bridge (Figure 3) based on normed distance (top) and statistical distance (bottom)

Figure 5 :

 5 Figure 5: Hierarchical clustering of modes obtained from clustering of 21 stabilization diagrams using the statistical distance with the same threshold D s lim = 16 as in the clustering of the diagrams (left), and the doubled threshold (right)

Figure 6 :

 6 Figure 6: Full stabilization diagram of Z24 Bridge

Figure 7 :

 7 Figure 7: Hierarchical clustering without post-processing of the full stabilization diagram of Z24 Bridge based on normed distance (top) and statistical distance (bottom)

Figure 8 :

 8 Figure 8: Hierarchical clustering after post-processing of the full stabilization diagram of Z24 Bridge based on normed distance (top) and statistical distance (bottom)

1

  Clustering Require: Frequencies f , damping ratios ζ, mode shapes ϕ, model orders n, maximum distance D lim , minimum number of modes per alignment nb stab Ensure: Clusters 1: Compute distance D ij between every pair of elements (Equation (7) or (9)) 2: Compute the clusters with complete linkage with respect to D lim 3: Sort the clusters by frequency and delete every cluster whose size is below nb stab Merge 4: Compare the neighboring clusters and merge the closest if their single linkage distance is below D lim Separation 5: for Every cluster do

	6:

ACKNOWLEDGMENTS

The support from the ANR "France Relance" program is gratefully acknowledged.