Guillaume Pichon

Alain Daidié

Éric Paroissien

Audrey Benaben

Clément Chirol

Experimental and Numerical Investigations of Mechanical Behaviour for Titanium Bolted Joints with Friction Shim

Keywords: Hole-to-hole, Fatigue, Static, Testing, Finite Element Analysis

The work presented here is an experimental and numerical study on a new technology called friction shims. Those are thin shim made of steel with a shape of washer. Their coating made of nickel and micro diamonds aim to increase the friction coefficient when located at the interface of two parts. The mechanical behaviour as well as the static and fatigue strengths of titanium bolted joints with this specific technology are presented in this study. Specimens are representative of hole-to-hole aeronautical joints, involving large clearance and misalignments. Specimen geometrical and material definition are presented as well as test set up. The experimental test results show limited impact on the static strength and a large increase in the fatigue strength. A numerical analysis based on the finite element method was developed and correlated with the experiments to provides complementary information and to understand the mechanical behaviour observed during the static tests and the gain in fatigue strength.

Introduction

One of the main concerns of aircraft manufacturers is to be integrated into industry 4.0, for better agility in both design and manufacturing. This transition is achievable by the re-engineering of processes. This paper is set in the framework of the "hole to hole" bolted assembly process introduced by Bloem [START_REF] Bloem | Developments in Hole-to-Hole Assembly[END_REF]. This process constitutes a great opportunity to revolutionize aeronautic assembly lines by having structural parts individually drilled before the assembly phase. However, the main consequence of this process is the potential for misalignments of bolt holes when parts are mated, making the bolt installation impossible. The solution proposed by Bloem [START_REF] Bloem | Developments in Hole-to-Hole Assembly[END_REF] is to slightly increase the diameter of bolt holes so that assembly remains feasible. The investigations on bolt hole clearance and error locations have been the subject of many publications, in particular concerning composite laminates [START_REF] Askri | Effect of Hole-location Error on the Strength of Fastened Multi-Material Joints[END_REF][START_REF] Lecomte | An analytical model for the prediction of load distribution in multi-bolt composite joints including hole-location errors[END_REF][START_REF] Liu | Probabilistic bolt load distribution analysis of composite single-lap multi-bolt joints considering random bolt-hole clearances and tightening torques[END_REF][START_REF] Mccarthy | Bolt-hole clearance effects and strength criteria in single-bolt, single-lap composite bolted joints[END_REF][START_REF] Lawlor | An experimental study of bolt-hole clearance effects in double-lap, multi-bolt composite joints[END_REF]. According to these latest publications, clearance and misalignments have a significant effect on static and fatigue strengths. Such defects also have an impact on the loading behaviour, in particular on the macro slipping, i.e. the slipping occurring at a joint interface when the shanks of the bolts and the holes of the plates are not in direct contact due to clearance. This macro slipping is an unwanted phenomenon. It is not present in current aeronautical joints due to fit clearance but it is inherent in hole-to-hole assembly. To delay this phenomenon, friction load at the joint interface is the best leverage to work on.

There are two ways to increase friction load at the interface. The first one is to increase the normal load, i.e. increase bolt loads. Several studies have investigated the effect of increasing bolt loads on the static and fatigue strengths [START_REF] Benhaddou | Effect of axial preload on durability of aerospace fastened joints[END_REF][START_REF] Mccarthy | An analytical model for the prediction of load distribution in highly torqued multi-bolt composite joints[END_REF] experimentally and numerically. Those works have shown a positive impact of increasing bolt load on fatigue life of bolted joints. However bolt load is limited to bolt axial strength in particular in the aeronautic field were bolts, mostly working in shear, have low axial strength for mass gain purpose.

The second way to increase friction load is to increase friction coefficient at joint interface. The proposed solution evaluated in this work is based on this method and uses friction shims industrialized by 3M. This technology is already used in the industry mainly to increase torque transfer in rotating systems (automotive, energy…) or to delete unwanted macro or micro sleeping This study aims to identify and quantify the effect of technology that can increase friction when located at the interface of a titanium bolted joint. It considers the mechanical behaviour as well as the static and fatigue strengths. To the author's knowledge it is the first time that friction shims are studied to evaluated their impact on static and fatigue strength of aeronautics bolted joints. A finite element (FE) analysis will be developed to increase the understanding of phenomena at play in this study and to evaluate the best way to model them.

2

Specimen preparation and testing procedure

Test specimens

For the experimental test campaign, the specimens were single-lap bolted joints made of titanium (TA6V) plates. Table 1 and Figure 1 give the dimensions of specimens, all dimensions being computed from the nominal bolt diameter (d=6.35mm). The real bolt diameter was 6.33mm. 𝑑 ̃ corresponded to the real hole diameter depending on the configuration. Fatigue and static specimens had different widths. For static loading specimens, a width of W=7d was used to insure bolt-bearing failure, while a width of W=4d was used for fatigue tests.

In the frame of hole-to-hole, three different configurations were used for the specimen geometry. Configurations A and B had aligned bolt holes with diameters of 6.38 and 6.88mm, respectively, while configuration C presented hole location errors in the loading direction, with opposite values of -0.5 and 0.5mm, and the diameter of the holes was 6.88mm. Table 2 presents the three different configurations and their associated clearances.

Friction washer

In this work, 3M industrial friction shims were used. Their geometry, positioning within the joint, and composition are given in Figure 2. A shim comprises a stack of three layers with a sheet of steel as the principal layer at the centre of the shim (0.15mm from micro measurement). On the exterior face, a thin layer of nickel is used to encapsulate diamonds. When the fasteners are installed, the diamonds create indentations in the specimen plates, thus creating higher friction grip. In this study grade 25 shims were used. The grade is based on the diamond average size, here 25µm. The shim has the shape of a washer with an inner diameter of 6.4mm and exterior diameter of 11mm. Corrosion risk are not a major concern at this step of the study. Furthermore titanium is subject to corrosion, but the oxide layer formed becomes very protective and waterproof and its thickness slowly increases.

Table 2. Configurations tested with clearance, misalignment and geometry

Testing setup

In this study, tensile tests to failure and fatigue tests were carried out to determine the joint behaviour in static and fatigue strength. Firstly, the tensile tests aimed at determining the static ultimate strength of the joint. It was evaluated by applying a constant displacement of 0.5mm/min at the tip of the specimen up to failure. The results were used to identify the stress range for the following fatigue tests and correlate the FE model in Section 4. Finally, the tensile tests were also used to assess the critical slipping load, characterized by a strong relative displacement between the two plates of the joint (so called macro slipping).

The fatigue tests were carried out with a loading ratio of R=σmin/σmax=0.1 and a frequency of 10 Hz. Each configuration was tested on seven specimens and the remote stresses applied were chosen to characterize the fatigue strength of the joints in the range of 10 4 to 10 6 cycles of fatigue life. The observed test results gave S-N curves and the Airbus Fatigue Indicator (AFI), corresponding to the stress level to reach a lifetime of 10 5 cycles.

Fig. 2. Friction shim positioning and composition from micro cut observations

3

Experimental test results and discussion

Static loading to failure

The static test results are gathered together in Figure 3. Those for configurations without friction shims are shown in dashed lines (A, B and C). The configurations with friction shims have solid lines (A_FS, B_FS and C_FS). The same nomenclature will be used for fatigue test results later in this paper. For reasons of confidentiality, the load has been normalized with the maximum ultimate load reached by configuration A_FS.

Joints without friction shim

At the beginning of the loading, all configurations follow the same path. In this phase, load is transmitted only by friction at the joint interface. Then, at a similar normalized load of 0.12, a change in curve slope can be observed for each specimen. The configurations A and C have a change of slope but the curve remains linear. However configuration B presents a nonlinear phase from 0.12 to 0.25 where the slope decreases then increases to finally become linear again with a slope similar to that of configuration A. This non-linear phase corresponds to macro-slipping. Theoretically, equation ( 1), can be used to compute the Coulomb's law friction coefficient, with FN the normal load, which is the bolt load, FT the tangential load and µ the friction coefficient. However, in this case, the macro slipping does not occur at a constant load level so it is hard to assess the correct tangential load. Assuming a friction coefficient of 0.4 for titanium/titanium contact [START_REF] Meier | In-process Measurement of the Coefficient of Friction on Titanium[END_REF] and knowing the normal load from previous instrumented installation tests, equation ( 1) gives a tangential normalized load value of 0.19. This value is in the middle of the macro-slipping zone and thus appears coherent.

FT =µ . FN (1)

It can be observed that this phase only exists for configuration B. For configurations A and C, the transition between the two linear phases is almost direct due to the geometry and clearance parameters of those configurations.

After the macro-slipping, the shank of fasteners and the holes come into contact and the load is then transmitted through the fasteners by bearing. This phase is linear while the deformation within the plate and the bolts remain elastic. This elastic deformation lasts until around 0.69 for configurations A and B. Finally, a nonlinear phase begins and lasts up to final failure of the specimensin this case by simultaneous shearing of the fasteners at a similar load of around 0.97.

For configuration C, the linear loading phase lasts until 0.46, when a nonlinear phase begins and continues until failure of the first fastener at 0.6 (fastener 2 in the scheme from Table 2.). At the instant of failure, the load drops from 0.6 to 0.29, then increases again until the final failure of the second bolt occurs at around 0.48.

Joints with friction shims

The global behaviour is very close to that without friction shims. The observed differences are described below.

The first linear phase, where load is transmitted by friction, reaches a much higher load: around 0.36 for configurations A and C, and 0.32 for configuration B.

After the slipping points for all the configurations, a small decrease of load, around 0.02, can be observed. This observation tends to indicate that the slipping behaviour of joints with friction shims has a static and a dynamic friction coefficient [START_REF] Berardo | An experimental-numerical study of the adhesive static and dynamic friction of micropatterned soft polymer surfaces[END_REF].

The macro-slipping phase is only visible on configuration B, similarly to a test without friction shim. However, it occurs at a quasi-constant load value of 0.30. Assuming the same bolt load, equation 1 gives values of 0.66 and 0.62 for the static and dynamic friction coefficients, respectively. Regarding the elastic deformation phase of the fasteners, the yield strength seems to be higher for joints with friction shims. However, the final failure reaches a similar load level of 1.00 for configurations A and B.

For configuration C, the maximal load reached before failure, occurs at the same moment, just before failure of the first bolt. However, the load level is higher for joints with friction shims, as it reaches 0.71.

Fatigue tests

Figure 4 presents the fatigue test results and the S-N curves extrapolated from raw data. For reasons of confidentiality, stress data have been normalized so that the Airbus Fatigue Indicator (AFI) of the configuration A is set to 1. The AFI corresponds to the stress value for a fatigue life of 10 5 cycles.

The results show two distinct phenomena. Firstly, at high stress levels, the more impacting parameter is the geometric configuration. As expected, configuration shows the highest fatigue life; for a normalized stress of its fatigue life is around 40000 cycles. Configuration B follows, with a very similar fatigue life. Finally configuration C shows the life, with 25000 cycles.

Secondly, at low levels of stress, the more impacting parameter is the presence or absence of the friction shim. The highest numbers of cycles reached before failure were observed at a representative stress value around 0.88 for the three configurations with friction shims. For the three configurations without friction shims, the highest number of cycles reached before failure were observed at representative stress around 0.7.

Finally, the three configurations with friction shims showed a 20% increase in their respective AFI compared with joints without friction shims. In particular, the AFI of configuration C_FS, despite the misalignment defects, was greater than the AFI from configuration A, representative of a conventional aeronautic joint.

Fatigue failure analysis

In this study, the design of specimens is such that failure occurs in the plates section close to the bolt holes. Due to the presence of two fasteners, the failure can occur on either the head part (hole 1), or the collar part (hole 2) as shown in Figure 5.

The failure location for each configuration is detailed in Table 3. It appears that the friction shim has no impact on the location on the failure. However, the geometric configuration tends to have an impact. For configuration A, the failure is equally distributed between the two sites; for configuration B, the failure occurs more often on the collar part; and for configuration C, the failure always occurs on the collar part. After observing the failure locations, fracture profiles were observed. In all cases, with or without friction shim, the fatigue failure started at the joint interface as shown in Figure 6. From those observations it can be noted that friction shims have no impact on the failure profile, despite the increased friction at the interface. Figure 6 also shows the gap between the two titanium plates when a friction shim is installed at the interface of the joint.

Finite Element analysis

To fully understands the observed experimental tests results, a three-dimensional (3D) FE model was set up with the commercial code Abaqus © . The geometrical, material, and stress parameters were respected according to their description in Section 2. Materials were defined with their elastic-plastic behaviour. The parameters remain confidential. The model has three calculation steps. During the first step, the bolts are loaded under tension, and both plates are clamped at their extremity with no displacement allowed. The second step is a step of relaxation. The last step is the loading step. Plate 1 (head side) is clamped at its extremity while a ramp of displacement is imposed at the extremity of plate 2 (collar side). The value of the displacements corresponds to the ones observed experimentally up to failure of the bolts, or to the failure of the first bolt for configurations C and C_FS. Contact is modelled with both a normal hard contact law and a tangential penalty law, for all the contact zone: bolt head/plate, collar/plate, bolt shank/plate, and plate/plate. A parametric study was run on the value of the penalty coefficient to match the experimental results. The final value of the penalty coefficient is 0.25 for configurations without friction shim. This value is lower than the one found experimentally, 0.4, but it is in agreement with previous results from McCarthy et al. [START_REF] Mccarthy | Experiences with Modeling Friction in Composite Bolted Joints[END_REF].

There were some convergence difficulties with configuration B. The large macro slipping due to the large clearance generated some computation errors. To stabilize the contact, the elastic slip of the tangential law was increased from 0.005 to 0.01.

To model the configurations with friction shim, the plate/plate contact region was decreased to match the surface area of the washers. In addition, the friction coefficient was strongly increased to 0.7, once again following a parametric study. The bolt loads remained unchanged as the same bolts were used in all of the configurations. Figure 9 presents the experimental and numerical load versus displacement curves. It can be observed that experimental and numerical results look very similar. The slight difference between A_FS and A_FS_NUM is due to the presence of clearance (0.06 mm) that was not taken into account in the numerical model. In the FE model, each contact load, normal or tangential, can be extracted. It was used to investigate the load distribution between fasteners and friction. In Figure 8, the evolution of each load: Bolt 1, Bolt 2 and Friction (with and without Friction Shim), is compared with the global load applied to the specimen for configurations C and C_FS. It can be observed that bolt 1 never transfers load. Static tests and numerical simulations have shown that failure of bolt 2 occurs before the contact of the first fastener, see Figure 7. These observations are coherent with the observed fatigue failure locations, as all failures were observed on hole number 2.

It can also be seen that the load carried by friction at the interface reaches a much higher maximum level for C_FS (0.28) than with C (0.1). Numerical results also show that the load transmitted by bolt 2 is higher by a value of 0.2 when the joint has no friction shims. Those observations are consistent with the experimental results.

Conclusions

In this study, an experimental campaign was conducted to investigate the effect of 3M friction shims at the interface of titanium bolted joints. The criteria studied were mechanical behaviour, quasi-static strength and fatigue life. A 3D numerical model by FE was developed to investigate the phenomena involved in greater depth.

Regarding the mechanical behaviour, the experimental tests showed that the macro slipping occurred at higher loads for joints with friction shims, and that the slipping behaviour was a combination of static and dynamic friction. Regarding static strength, friction shims have no impact on joints with aligned holes and have a low impact on joints with misaligned holes. The numerical model showed that, as the load increases in the bolts, the friction load decreases, explaining the observed results. In particular, for the configurations with misaligned holes, the heterogeneous distribution of load among the fasteners generated damage initiation in one of the bolts. The remaining one was not damaged, thus keeping the friction load effective.

Finally, regarding fatigue life, friction shims led to an increase of 20% on the Airbus Fatigue Indicator, even for joints with geometrical defects. The numerical model showed that, at the same loading level, the joints with friction shims were transmitting less load through the bearing of the bolts, thus reducing the stress in the plates.

The next steps on this project are to run similar experiments with different plate material in particular aluminium and composite material widely use in the aeronautic industry. Corrosions risks will also be evaluated more intensively reading aluminium as this material is more sensitive to corrosion than titanium.

Fig. 3 .

 3 Fig. 3. Normalized load vs. global displacement during static tests and zoom on the early loading phases

Fig. 4 .Fig. 5 .Fig. 6 .

 456 Fig. 4. Fatigue test results and associated S-N curves

Fig. 7 . 1 Fig. 8 .

 718 Fig. 7. Cut view of configuration C, under quasi-static loading at failure of bolt 2, Von Mises equivalent stress with deformation coefficient of 1

  was possible to obtain similar results numerically by strongly increasing the friction coefficient and reducing the elastic slip.

Fig. 9 .

 9 Fig. 9. Comparison between experimental and numerical results for configurations without friction shim on top and with friction shim on the bottom

Table 1 .

 1 Specimen geometry calculation and values

		d	𝒅 ̃	a	L	e	p		W	t
	Calculation	-	-		18d+a 2d	4d	7d (static)	4d (fatigue)	d/2
	Value (mm) 6.35	6.38 (A), 6.88 (B, C)	70 184.3 12.7 25.4 44.45	25.4	3.175

Fig. 1. Specimen geometry definition, numerical model and in reality


Raul Carreira
email: raul.carreira-rufato@espci.fr

Y Diouane

Joël Henry

Richard Ahlfeld
email: richard@monolithai.com

J Morlier
email: joseph.morlier@isae-supaero.fr

A mixed-categorical data-driven approach for prediction and optimization of hybrid discontinuous composites performance

Surrogate models are an essential engineering tool and their popularity has increased recently due to the high computational cost of evaluating real-world simulations. However, most of these functions are described by mixed variables (continuous and categorical), which makes it harder to create accurate interpolation functions. This work builds a surrogate model from a given mixed data set, in order to quickly and accurately calculate the mechanical performance of hybrid discontinuous composites. Then, in order to find the optimal hybridization, three different approaches are performed: mono-objective, targeted and multi-objective. Starting from a virtual database provided by the industrial partner, the mixed categorical optimization process is performed by coupling a multi-armed bandit strategy with a continuous Bayesian optimization solver. The efficiency of the proposed approach is tested and two main results are achieved. The obtained surrogate models are shown to be sufficiently accurate, having an 2 score grater than 90% in average. Our proposed optimization process is also able to identify correctly the optimal fibres with respect to the desirable targets.

I. Introduction

T , the majority of real-world simulations depends on specific categories and continuous design variables [1].

For example, a possible application is to predict the mechanical performance of hybrid discontinuous composite materials (i.e., materials manufactured by combining two or more different types of fibres). Most of these real-world black box systems can be extremely expensive to evaluate. For that reason, surrogate models are often used to accelerate the design process. Therefore, it is necessary, for engineering tasks, to build an accurate surrogate model which deals with both categorical and continuous variables and facilitates the decision-making in engineering design process.

It is commonly known to use categorical variables as outputs of machine learning models classification-based (i.e., classes). However, the work presented in this paper deals with the categorical variables as inputs to the prediction model, which often requires different methodologies. Various reviews of existing mixed continuous and categorical surrogate modeling techniques exists [START_REF] Herrera | Metamodel-assisted optimization based on multiple kernel regression for mixed variables[END_REF][START_REF] Halstrup | Black-box optimization of mixed discrete-continuous optimization problems[END_REF][START_REF] Roustant | Group kernels for Gaussian process metamodels with categorical inputs[END_REF]. One approach [START_REF] Herrera | Metamodel-assisted optimization based on multiple kernel regression for mixed variables[END_REF] consists in converting the categorical variables into continuous (e.g., dummy coding) and including them in standard continuous surrogate models. Another well-known approach is based on the construction of a concept of distance between mixed variables and then, on the construction of a surrogate model through a Gaussian process [START_REF] Halstrup | Black-box optimization of mixed discrete-continuous optimization problems[END_REF]. A recent work [START_REF] Pelamatti | Overview and Comparison of Gaussian Process-Based Surrogate Models for Mixed Continuous and Discrete Variables: Application on Aerospace Design Problems[END_REF] describes and compares the adaptation of surrogate modeling techniques based on Gaussian processes for mixed continuous-categorical variables. In particular, the authors in [START_REF] Roustant | Group kernels for Gaussian process metamodels with categorical inputs[END_REF] investigate the well-known kernels group, and study the generalized compound symmetry covariance matrices, that are parsimonious parameterizations of the covariance matrix, written in block form with a fixed covariance between pairs of blocks and within blocks.

The use of data-driven approach in surrogate models is a technique used when it is not possible (or expensive) to have access to the simulation code. Recently, the authors in [START_REF] Finley | Data-driven intelligent optimisation of discontinuous composites[END_REF] used this concept of a data-driven intelligent routine to find the best aligned discontinuous composites. The obtained results showed the potential of the proposed approach in a single and multi-objective optimization contexts. The incorporation of categorical inputs, in addition of the continuous ones, can be done using a random forest methodology [START_REF] Hutter | Sequential model-based optimization for general algorithm configuration[END_REF]. However, such approach has some limitations in performing extrapolation to unknown zones of the design space. Recently, the authors in [START_REF] Nguyen | Bayesian Optimization for Categorical and Category-Specific Continuous Inputs[END_REF] proposed a method that uses a multi-armed bandit strategy to handle categorical inputs. In this proposed approach, each category corresponds to an arm with its reward distribution centered around the optimum of an objective function with continuous inputs. The obtained results [START_REF] Nguyen | Bayesian Optimization for Categorical and Category-Specific Continuous Inputs[END_REF] were very promising and show that the proposed method is able to explore well the categorical design space with a minimal computational cost. For that reason, this research chose to work with this method to perform the optimization process.

The use of machine learning (ML) algorithms in mechanical engineering has been studied a lot recently, Geetha et al. [START_REF] Geetha | Overview of machine learning and its adaptability in mechanical engineering[END_REF] presents an overview of some applications. Additionally, it has been also used to help in material design decision-making process [START_REF] Huang | Artificial Intelligence in Materials Modeling and Design[END_REF]. Interestingly, the authors in [START_REF] Bessa | A framework for data-driven analysis of materials under uncertainty: Countering the curse of dimensionality[END_REF] build a new data-driven computational framework to help design new materials and structures, building their own training set to represent the macroscopic response. Using supervised learning technics, Nyshadham et al. [START_REF] Nyshadham | Machine-learned multi-system surrogate models for materials prediction[END_REF] demonstrates the utility of surrogate models to interpolate the characteristics of the atomic structure for different alloys, showing that the surrogate models can be used to select materials and their most desirable properties in a wide range of possibilities. Surrogate models are also used in [START_REF] Hou | A novel development of bi-level reduced surrogate model to predict ductile fracture behaviors[END_REF] to identify the material parameters of ductile fracture criterion. Recent works, [START_REF] Chang | Prediction of the Composition and Hardness of High-Entropy Alloys by Machine Learning[END_REF][START_REF] Liu | Machine learning-assisted modeling of composite materials and structures: a review[END_REF] present how ML algorithms can facilitate the prediction of important material properties.

In the context of composite materials, ML approaches have been widely used [START_REF] Bessa | A framework for data-driven analysis of materials under uncertainty: Countering the curse of dimensionality[END_REF][START_REF] Nyshadham | Machine-learned multi-system surrogate models for materials prediction[END_REF][START_REF] Pathan | Predictions of the mechanical properties of unidirectional fibre composites by supervised machine learning[END_REF][START_REF] Pattnaik | Role of machine learning in the field of Fiber reinforced polymer composites: A preliminary discussion[END_REF][START_REF] Yang | Prediction of composite microstructure stress-strain curves using convolutional neural networks[END_REF]. Chen et al. [START_REF] Chen | Machine learning for composite materials[END_REF] discuss and review recent progresses in ML applied to composite materials design, specifically the application of those models in the prediction of the materials' mechanical properties. Additionally, the authors in [START_REF] Pathan | Predictions of the mechanical properties of unidirectional fibre composites by supervised machine learning[END_REF] train the gradient-boosted tree regression model to predict the macroscopic stiffness and yield strength of unidirectional composites. Their method is interesting as it analyses the material's micro-structure image without performing physically-based simulations. These works were able to demonstrate the excellent results in the use of machine learning models in predicting the characteristics of composite materials in order to help in engineering design. However, limitations still lie in the fact that some characteristics of hybrid composite materials can be described by categorical variables, which are usually handled by some coding method. These methods can significantly increase computational time as well as limit correlation and neighborhood characteristic analyzes between categorical and continuous variables. Our work comes precisely to contribute in this aspect. In this paper, using a mixed categorical-continuous database, we construct mixed surrogate models and then use optimization approaches to make final recommendations to help the complex engineering design process. To the best of our knowledge, this paper presents the first mixed categorical-continuous data-driven approach to predict and estimate optimal mechanical performance for hybrid discontinuous composites.

Our proposed method is used and validated in the determination of the characteristics of hybrid discontinuous composites. In fact, due to the wide variety of possible combinations of different types of fibres, it is becoming increasingly not viable and very expensive to carry out stress and deformation tests for this type of material. Thus, the proposed method for optimization and construction of accurate surrogate models is of high interest. This work, in particular, helps to optimize the scalar characteristics of the stress-strain curves for hybrid discontinuous composite materials (combination of two different fibre types). The main characteristics of this type of material can be seen in Fig. 1. The applications of composite materials are increasing in several engineering areas. Composite structures are usually chosen according to specific design goals, due to the individual constituent fibers and their dimensions and directions. Therefore, to obtain the best results, optimization models are used to find the best solutions that satisfy a given set of design constraints. [START_REF] Maalawi | Introductory Chapter: An Introduction to the Optimization of Composite Structures[END_REF] presents a review on the optimum design of composite structures and the relevant optimization techniques. Various optimization methods are classified according to the optimization process, but the most common approaches are conjugate directions and conjugate gradient methods. Moreover, genetic algorithms, which are commonly used methods, are part of a group of methods that do not need information about the function's gradient to find the global optimum and are generally used in optimization of composites structures [START_REF] Narayana Naik | Design optimization of composites using genetic algorithms and failure mechanism based failure criterion[END_REF]. For [START_REF] Maalawi | Introductory Chapter: An Introduction to the Optimization of Composite Structures[END_REF], genetic algorithms has been the most efficient method for obtaining the global optimum design of composites. In this paper, we used another stochastic method for the optimization that was implemented by [START_REF] Nguyen | Bayesian Optimization for Categorical and Category-Specific Continuous Inputs[END_REF].

For the purpose of this work, we will first build a surrogate model (on the mixed categorical design space) and then optimize the output of the obtained surrogate. The main goal is to find the best fibre combinations according to the desired scalar parameters.

This paper is organized as follows. Section II presents the theoretical background required by our methodology to solve the regarded problem. The application problem, with its formulation and the construction of the database, is detailed in Section III. In Section IV, we discuss the obtained results. Final conclusions are drawn in Section V.

II. Theory and methodology elements

Much of today's physical phenomena are given as expensive-to-evaluate black-box simulations. In this context, building explicit inexpensive surrogate models can be very useful to accelerate the physical modeling process. Surrogate models ought to build approximate (often polynomial) functions that interpolate the black-box simulation function and extrapolate the behaviour to unknown zones of the design space. Gaussian process are an efficient tool to build cheap and accurate surrogate models.

A. Gaussian Process

Gaussian Process with real inputs

Gaussian Process (GP) [START_REF] Abrahamsen | A review of Gaussian random fields and correlation functions[END_REF][START_REF] Cressie | Fixed rank Kriging for very large spatial data sets[END_REF][START_REF] Hensman | Gaussian processes for big data[END_REF][START_REF] Matheron | La Théorie des Variables Régionalisées et ses Applications[END_REF], also known as Kriging models, return a meta-model function of the form ( ) ∼ N ( ( ), ( ) 2 ) for any input ∈ R , where represents the mean (deterministic) and is the standard deviation that depicts the GP uncertainty of each sample on the entire domain.

A detailed description of GPs with real inputs can be as follows. Let : R ↦ → R be a scalar function for which a GP is built using a DoE of points D = { , } =1,..., where ∈ R and = ( ) ∈ R. The GP model, related to using sample points, is a family of functions defined by a mean function and a standard deviation . Namely, at each point , the GP of is defined with a multivariate Gaussian distribution N ( ), ( ) 2 . In the context of this paper, the mean and the standard deviation are given by Equation 1

( ) = ( ) ˆ + ( ) -1 -ˆ ( ) = ˆ 1 -( ) -1 ( ), (1)
where : R ↦ → R is a vector of polynomial basis functions, = [ ( 1 ), . . . , ( )] ∈ R × is the matrix of predictors at the sampling points in D and = [ 1 , . . . , ] ∈ R is a vector of outputs of . For a given kernel function : R ×2 ↦ → R, we set ( ) = [ ( 1 , ), . . . , ( , )] ∈ R as the correlation vector between the point and the sample points of D. is the correlation matrix = [ ( , )] , =1,..., ∈ R × computed over all the sample points of D. Finally, the regression coefficients are defined by Equation 2

ˆ = -1 -1 -1 ∈ R , (2)
and the the standard deviation factor ˆ of the GP is given by Equation 3

ˆ = 1 -ˆ -1 -ˆ . (3)
We note that the means and are computed thanks to a correlation function chosen by the user (see the next subsection). In fact, the definition of the kernel function may depend on a set of hyper-parameters that are in general estimated by maximizing a likelihood function.

Gaussian Process for mixed inputs

An efficient extension of GP to handle mixed (categorical and continuous) input variables was proposed in [START_REF] Pelamatti | Overview and Comparison of Gaussian Process-Based Surrogate Models for Mixed Continuous and Discrete Variables: Application on Aerospace Design Problems[END_REF][START_REF] Roustant | Group kernels for Gaussian process metamodels with categorical inputs[END_REF]. A detailed description of the proposed methodology is as follows. Consider : R × F ↦ → R be a scalar function for which a GP is built using a DoE of points of the form D = { , , } =1,..., where ∈ R , ∈ F and = ( , ) ∈ R. The space F is a set of categorical variables 1 , ..., with 1 , ..., levels, respectively. The number of levels results in a total number of categorical combinations of =1 . The kernels on D can be obtained by combining kernels for continuous and categorical variables. Let be the correlation kernel of continuous variable inputs and be the one for the categorical inputs. Different choices for the global correlation matrix , with respect to mixed inputs, are possible. For example, can be set to .

,

+ or to (1 + ) (1 +
). In this paper, following the suggestions in [START_REF] Roustant | Group kernels for Gaussian process metamodels with categorical inputs[END_REF], we set, for all , = 1 . . . , , ( ,

) = ( , ) ( , ), (4)
where = ( , ) for all = 1, . . . , . The correlation matrix for the mixed inputs is then set to be

= [ ( , )] , =1,..., . (5)
For the continuous inputs, we will consider exponential kernel adapted to Gower distance correlation kernel, i.e.,for all , = 1 . . . , ,

( , ) = 2 exp - =1 | -| Δ , (6)
where, for all = 1, . . . , , > 0 and 1 ≤ ≤ 2, are hyper-parameters, 2 is the variance associated with the continuous kernel. For all = 1, . . . , , Δ is an estimation of the range, i.e., Δ = ( + ) max =1,...

min =1,...

. For the categorical inputs, the distance is measured based on a similarity score combined with the p-exponential kernel [START_REF] Halstrup | Black-box optimization of mixed discrete-continuous optimization problems[END_REF][START_REF] Gower | Euclidean distance geometry[END_REF], i.e., for all , = 1 . . . , ,

( , ) = 2 exp - =1 ˆ ( , ) + ˆ , (7)
where, for all = 1, . . . , , ˆ > 0 and 1 ≤ ˆ ≤ 2, are hyper-parameters, 2 is the variance associated to the categorical kernel. is a score function given by

( , ) = 0, if = 1, if ≠ . ( 8
)
Similarly to the continuous case, the correlation vector between a given point = ( , ) and the sample points of D = { , } =1,..., is given by

( ) = [ ( 1 , ), . . . , ( , )] . (9)
For the sake of simplicity, in this paper, we consider to use a constant trend (zero order polynomial) for the polynomial basis functions is ( ) = 1 ∈ R. Hence, = 1 ∈ R a vector of size where all the components are set to 1. For a given ∈ R × F , the GP meta-model function

( ) ∼ N ( ( ), ( )
2
) where and are computed using the same formulas used for continuous inputs (given in Section II.A.1) where , ( ), ( ) and are replaced with , ( ), ( ) and , respectively. We note that, as in the continuous case, the GP models for mixed inputs will require the estimation of 2( + ) + 2 hyper-parameters (i.e., , , , ˆ , and ˆ ) which is done by maximizing a likelihood function. In this paper, we considered that = ˆ = 1 and = = 0 , so only the hyper-parameters ( 0 , , ˆ ) are estimated. Now that our surrogate model has been built using a mixed training set D . In the next section, we will detail our proposed optimization approach to perform predictions of unknown combinations of continuous-categorical inputs, possibly not present in our initial training set D .

B. Bayesian optimization

Bayesian optimization [START_REF] Frazier | A Tutorial on Bayesian Optimization[END_REF] (BO) is an extremely useful tool to handle black-box and expensive-to-evaluate optimization problems. BO targets to find the global optimum * (a maximizer or a minimizer) of an unknown objective function . Namely,

* = arg max ∈Ω ( ), (10)
where Ω is the design space of interest. We note that, for simplicity reasons, we consider only a maximization formulation. A minimization problem with the objective functioncan be seen also as a maximization problem.

Bayesian optimization for continuous inputs

The BO method with real inputs [START_REF] Jones | Efficient global optimization of expensive black-box functions[END_REF] can be described as follows. Let : R ↦ → R be a scalar function, known as the objective function. Let D = { , } =1,..., be a design of experiments (DoE) of size such that ∈ R and = ( ) ∈ R for all = 1, . . . , . Let and be the mean and the standard deviation of the GP associated with the DoE D (as in Section II.A.1).

BO is then a sequential search process in which, at iteration , one enrich the DoE using a new point ( +1 , +1 = ( +1 )). The enrichment process is led by an acquisition function ( ) modeling the trade-off between exploration of new areas in the design space (i.e., areas with high value of ) and exploitation (i.e., minimization of ). Some well-known acquisition functions can be find here [START_REF] Jones | Efficient global optimization of expensive black-box functions[END_REF][START_REF] Kushner | A new method of locating the maximum point of an arbitrary multipeak curve in the presence of noise[END_REF]. After max objective function evaluations, the algorithm return a final recommendation max . A brief description of the BO process is given in Algorithm 1.

Algorithm 1: BO for continuous inputs.

Data: : The initial DoE D 0 . A maximum number of iterations max ; for = 0, . . . , max -1 do

Step 1: Build the surrogate models using GPs using the DoE D ;

Step 2: Find +1 a solution of the enrichment maximization sub-problem;

Step 3: Evaluate the objective function at +1 = ( +1 ) and update the DoE: D +1 ; end Return: The best point found in the final DoE D max ;

Bayesian optimization for mixed inputs

To handle the mixed categorical-continuous inputs, different BO based approaches can be used. One possible approach is by applying the framework given in Algorithm 1 using GP models for mixed inputs (as given in Section II.A.2). This would mean that one needs to solve the enrichment maximization problem (Step 2 of Algorithm 1) in the mixed input space R × F . Due to the combinatorial nature of F , exploring the full mixed space can be difficult, or even out of reach for large values of and .

In this work, we consider to use a different approach, where the enrichment maximization problem (Step 2 of Algorithm 1) will be performed only over the continuous design variables [START_REF] Nguyen | Bayesian Optimization for Categorical and Category-Specific Continuous Inputs[END_REF]. This will allow to avoid the curse of dimensionality (related to the categorical inputs) and thus the maximization problem in Step 3 will be solved with a minimal computational cost.

The regarded approach formulates the problem as a multi-armed bandit (MAB) problem [START_REF] Bubeck | Regret Analysis of Stochastic and Nonstochastic Multi-armed Bandit Problems[END_REF][START_REF] Slivkins | Introduction to Multi-Armed Bandits[END_REF], where each category corresponds to an arm with its reward distribution centred around the optimum of the objective function over the continuous input space. Then, it identifies the best arm and maximizes the corresponding continuous function simultaneously. This is possible by considering each arm of the method as a combination of different categories. Finally, the approach uses a Thompson sampling scheme [START_REF] Thompson | On the likelihood that one unknown probability exceeds another in view of the evidence of two samples[END_REF] to connect the MAB approach (used to explore the categorical space) and the Bayesian optimization approach (used to handle the continuous inputs) in a unified framework. A detailed description of the optimization is given in [START_REF] Nguyen | Bayesian Optimization for Categorical and Category-Specific Continuous Inputs[END_REF].

III. Industrial test case and overall approach

The methodology presented in Section II will be applied to a real industrial context where, using a virtual database D , the main goal is to build a prediction approach for hybrid discontinuous composite material. The database used in our study was generated by the "Virtual Testing Framework" [START_REF] Finley | The influence of variability and defects on the structural performance of complex composite microstructures[END_REF][START_REF] Henry | Semi-analytical simulation of aligned discontinuous composites[END_REF] and was created to optimize the mechanical properties of hybrid discontinuous composites [START_REF] Finley | Data-driven intelligent optimisation of discontinuous composites[END_REF]. The regarded framework considers a hierarchical approach where the interactions between fibres at micro scale are first considered and then used at macro level (see Fig. 2).It is out of scope of this work to explain deeply the development of this database, additional information can be obtained regarding [START_REF] Finley | The influence of variability and defects on the structural performance of complex composite microstructures[END_REF][START_REF] Henry | Semi-analytical simulation of aligned discontinuous composites[END_REF].

In this context, the mechanical response of interactions between fibres of different types (diameter, stiffness, strength...) was calculated using a mathematical shear lag model developed in [START_REF] Henry | Semi-analytical simulation of aligned discontinuous composites[END_REF] and extended to hybrid microstructures in [START_REF] Henry | Modelling hybrid effects on the stiffness of aligned discontinuous composites with hybrid fibre-types[END_REF][START_REF] Henry | Virtual testing framework for hybrid aligned discontinuous composites[END_REF] (see Fig. 2a). The stress-strain curves of each fibre is then computed by averaging the stress-strain curves of interactions and adding fibre breakage (Figs. 2b and2c). In Fig. 2c, one can see that fibre breaks will occur as the stress increases. The overall stress-strain curve of the composite is then computed (Fig. 1). Finally, the final failure is calculated from a failure criterion [START_REF] Henry | Semi-analytical simulation of aligned discontinuous composites[END_REF] based on Dugdale's non-linear fracture mechanics criterion [START_REF] Dugdale | Yielding of steel sheets containing slits[END_REF]. This Framework enables to accurately estimate the stress-strain curve of a hybrid aligned discontinuous composite of any size.

Our test database consists of approximately 10 4 points where each point of the database is described by seven input variables (2 categories and 5 continuous values) and five scalar output variables. The five continuous are described in Table 1, see also Table 2 for a sample of 5 points extracted from the tested database.

Continuous inputs

: Length of the fibres ( ) : Proportion of carbon fibres in the composite S : Shear strength of the matrix (MPa) G : Tangent shear stiffness (MPa) F : Fracture toughness of the matrix ( / 2 ) Table 1 A detailed description of the continuous inputs.

The categorical input variables are sets of different types of glass fibres and carbon fibres. We note that in this work, the words "glass" and "carbon" actually mean "low stiffness" and "high stiffness". So the same fibre can be considered 3 A detailed description of the categorical sets: carbon and glass fibres.

The aim of our study is twofold. First, we will build surrogate models for this hybrid discontinuous composite materials (using a mixed database D ), as explained in Section II.B.2. Second, we will make a final recommendation by performing 3 different mixed categorical-continuous BO optimization approaches: mono-objective, multi-objective and inverse optimization.

The overall approach proposed in this paper is outlined in the flowchart given by Fig. 3. Using a mixed database D , we first build a set of five surrogate models using GP. For all = 1, . . . , 5, let : R × F → R be the GP model associated with the th output. The categorical set F combines all the possible arrangements of carbon and glass fibres in our database (namely, = 15 × 16 = 240). We note that, in our database, we dispose of 5 continuous inputs, hence = 5. Additionally, for practical reasons, the search space of continuous inputs R will restricted to the set [ , ] := { ∈ R , ≤ ≤ } where ∈ R and ∈ R are respectively the lower and the upper bounds. In what comes next, we will use Ω to denote the mixed targeted search space, i.e.,

Ω := = ( , ) ∈ R × F | ∈ [ , ] . (11)
Once the surrogate models are built, as outlined in Fig. 3, one can start the optimization process following 3 different approaches: Approach I: a mono-objective approach. In this approach, we will try to maximize one desired output (e.g., the initial stiffness). Namely, we select an output of interest and perform optimization of over Ω, i.e., max ∈Ω ( ) := ( ).

Approach II: an targeted approach. In this approach, we will try to find the set of inputs that allow to fit the best a set of targeted (or observed) output data. Namely, given a set of targeted output data ( * 1 , . . . , * ) ∈ R where, for all ∈ {1, . . . , }, ∈ {1, . . . , 5} (i.e., a selected subset of the output scalars, see Table 2). Our goal is to find the input variables that minimize the misfit function. In this case, we will solve the following normalized maximization problem:

max ∈Ω ( ) := - =1 ( ) - * 2 1 + ( * ) 2 . ( 12
)
Approach III: a multi-objective approach. This approach is based on a concurrent optimization process which allows experts to incorporate multiple conflicting objectives and to specify the trade-offs between them. Instead of obtaining one single optimal solution (as given by Approach I), a set of non-dominated solutions (known as the Pareto front) can be proposed. Using a linear scalarization approach it is possible to approximate the Pareto front. Namely, for each :

= { } 1≤ ≤ such that { =1 = 1 | ∈ [0, 1]},
where ( 1 , . . . , ) are the selected output surrogate models. It is important to emphasize that the choice of the optimization method depends on the type of problem encountered. Therefore, the three different methods presented in Fig. 3 are not confronted, but must be chosen initially depending on the desired optimization type.

IV. Results and discussion

A. Surrogate model results

Using our database, we built mixed surrogate models for the outputs. The surrogate modeling process, as described in Section II.B.2, was implemented within the surrogate modelling toolbox (SMT) [START_REF] Bouhlel | A Python surrogate modeling framework with derivatives[END_REF]. SMT is an open-source Python package consisting of libraries of surrogate modelling methods with only continuous inputs.

In order to compare the efficiency and robustness of our method, a comparison was made with 2 other existing models applied to our database. The first model used in the comparison was the one described in [START_REF] Finley | Data-driven intelligent optimisation of discontinuous composites[END_REF], in which a technique for encoding categorical variables called one-hot-encoding was used. This way of coding the categorical variables significantly increases the input dimension by one extra variable per category, because all the extra variables are set to zero except the required category which is set to one. After the codification process, this approach treats extra variables as continuous between 0 and 1. Then a Gaussian regression model of continuous variables like the one described in Section II.A.1, using the Marten 5/2 kernel, was created. The second method used for comparison, consisted of a previous codification of the categorical variables in integer values from 1 to the number of categories, thus excluding the notion of distance between categorical variables. Once the variables were coded, a random forest model was used to generate the surrogate model.

Models

Pearson coefficient ( 2 score) Due to the large size of our database, we used only 1000 points (selected randomly) to train and test the surrogate model faster. It is noteworthy that due to the increased dimension of the model when using one-hot-encoding, the computational training time is much higher than the other two types of model. The presented results are displayed using the average of 10 different runs. For each run, our data set was the split to two parts, 75% for training and 25% for testing the model. The obtained results for the five outputs for three different models are given in Fig. 4. For each scalar output, the Pearson coefficient 2 is also calculated, The 2 score measures the linear correlation between the predicted and the true values (Tab. 4 shows the 2 score for the three different models). It is still possible to realize that using only 1000, our model was the one that got the best performance in relation to the 2 . However, one can see, from the obtained values of 2 of our model, that the prediction accuracy (although high) was not the same for all the scalar outputs.

This can be easily explained by the stochastic nature of each output. The strain outputs (top row in Fig. 4) are much harder to predict because their variability is very high as they are a "weakest-link-based" variable: an apparently similar composite could fail at very different strains. Strength predictions (second row) are also affected by variability in the strength, although this effect is less strong. Finally, the stiffness (third row) is an "average-based" variable, and is therefore less prone to vary than failure variables. This makes it easier to predict, and explains that the accuracy of stiffness prediction is the highest. In any case, the displayed 2 scores of our model were all above 90%, which indicates that the constructed models are sufficiently accurate.

A sensitivity analysis was performed to evaluate the impact of the inputs on the accuracy of the constructed surrogate models ( 1 , 2 , ... , 5 ). The obtained 2 scores are presented in Table 5. Our strategy here consisted of discarding a subset of inputs from the database and training the surrogate model for each situation, without such inputs. For each output, the obtained Pearson coefficient will be compared with its original value when all inputs were used.

Table 5 shows clearly that some inputs do not have a high impact on the accuracy of the prediction result. In this case, such inputs can be neglected and this, in particular, will allow later to reduce the computational cost of the optimization process (see the next subsection). In fact, by suppressing the three least influential inputs (namely, the input variables S, G and F), the Pearson coefficient for the obtained surrogate models are at worst 6.6% lower than the models built using all the inputs. As a consequence, in what comes next, we will restrict the definition of the mixed search space Ω to be of the form: where ∈ R 2 and ∈ R 2 .

Ω := = ( , ) ∈ R 2 × F 240 | ∈ [ , ] , (14)

B. Optimization results

Using the trained surrogate models, we will now try to perform the optimization process following the three approaches as given in Section III. The obtained recommendations will be compared to desirable targets specified by the industrial partner "Monolith AI". Table 6 details the desirable targets and the associated optimization strategies.

For the BO solver we used the Python library Bandit-BO * which implements the approach described in Section II.B.2. All the parameters were kept unchanged, except the maximum number of iterations (set to 100), the lower bound (set to [0, 515] ) and the upper bound (set to [START_REF] Pelamatti | Overview and Comparison of Gaussian Process-Based Surrogate Models for Mixed Continuous and Discrete Variables: Application on Aerospace Design Problems[END_REF]12000] ). Due to the stochastic natures of the Bayesian optimization solver, we performed 10 different runs and display the results in average. more advanced mixed-categorical GP models can lead to significant improvement in our results, for instance by using continuous relaxation [START_REF] Saves | Multidisciplinary design optimization with mixed categorical variables for aircraft design[END_REF]. These ideas as well as the inclusion of extensive numerical tests will be addressed in a forthcoming work.

Fig. 1 A

 1 Fig. 1 A typical stress-strain curve for a hybrid discontinuous composite material.

  central fibre (red type) interacts with its four neighbours. (c) Interaction between two types of fibres.

Fig. 2

 2 Fig. 2 Overview of the different scales present in the "Virtual Testing Framework" for aligned hybrid discontinuous composites [34].

Fig. 3

 3 Fig. 3 An overview of the data-driven optimization process as proposed in this work.

  we solve the following parametrized (by ) mono-

Fig. 4

 4 Fig. 4 Obtained results for the five scalar outputs. The black line represents the perfect scenario where the predicted values are equal to the true ones.

Table 2 A sample of 5 points extracted from our experimental database. as

 2 one or the other depending on what is the other fibre. The set of carbon fibres is composed of 16 different levels while the set of glass fibres is composed of 15 different ones. A detailed description of the levels is given in Table3.

	Carbon Fibres
	XN-90	P120J	T1000GB	C124
	T800H	M60JB	C320	M40B
	P75S	K13D	T300	XN-05
	FliteStrand SZT HTA5131 GF	C100
	Glass Fibres
	GF	XN-05	FliteStrand SZT C124
	T300	T800H	C320	P75S
	C100	XN-90	HTA5131	T1000GB
	P120J	M40B	M60JB
	Table

Table 4 A comparison between three different models: Mixed GP (our model), One Hot Marten ([5]) and a Random Forest model. In bold the scores below 90%. (average of 10 runs)

 4

		1	2	3	4	5
	Mixed GP	0.994 0.925 0.935 0.959 0.972
	One Hot Marten 0.992 0.726 0.899 0.900 0.911
	Random Forest	0.900 0.831 0.935 0.805 0.860

* https://github.com/nphdang/Bandit-BO

Acknowledgments

We would like to thanks the "Fondation ISAE-SUPAERO" for their financial support.

Inputs

Pearson coefficient ( 2 score) For the problems I(a), I(b) and I(c), only 100 iterations for the BO solver were used to reach convergence. For the problems II(a), II(b) and II(c), the optimization process was harder, a specific DoE was used and 500 iterations for the BO solver were used to obtain the desired convergence. The obtained convergence plots are depicted in Figure 5. Table [START_REF] Nguyen | Bayesian Optimization for Categorical and Category-Specific Continuous Inputs[END_REF] shows the obtained optimization results associated with the same optimization process. 6 The desirable targets of the industrial partner "Monolith AI" and the associated optimization strategies.

The names of the problem depend on the optimization approach we consider to solve it (I, II, or III).

that the Tab. [START_REF] Nguyen | Bayesian Optimization for Categorical and Category-Specific Continuous Inputs[END_REF] does not present the results for the multi-objective optimization (Problem III), since the best choice depends on the criterion adopted and the variable to which it seeks greater importance. The results for this case are presented by the Pareto front in of Fig. 6. Clearly, for all the runs, the categorical choices were almost all compatible, the expected values appeared with higher occurrence. It can be highlighted that in case II(a), K13D and XN-90 have similar behaviours and generally, when one fibre type provides good results, the other one should be good as well. Moreover, the "expected value" column was not necessarily listing all the values expected, but at least some that were expected. Then if the optimization finds other values, it should not be a problem, it just means that in that specific case, the engineer will have the choice between two different options, because they would have a similar behavior. Therefore, in case II(a) the choice could be done by the price of the fibre, choosing the cheapest one to reduce cost. The scalar values were slightly different following two scenarios. In the first one, the values were very close to the predicted range. We believe that such values can be exactly reached using a higher number of iterations within the BO solver. The second scenario is related to outliers where possibly errors are propagated during the construction of the surrogate GP model.

In the second part of the analysis, we tried to solve problem III(a) and III(b) (see Table 6). The associate Pareto Fronts are depicted in Fig. 6. For each point of the Pareto Front we run 10 scalar optimization problems. The mean and the standard deviation are presented for each point.

Due to the categorical nature of the input variables, the Pareto fronts present a clear gap. This is in particular due to the absence of intermediate mixed variables in the Pareto-optimal region. Regarding the categorical choices, clearly the two obtained Pareto front are each clustered around three (Carbon, Glass) fibres. For the different runs, the Pareto front associated with the problem III(a) (where we maximize the initial stiffness and the pseudo ductile strain, see Fig 6a Fig. 6 The obtained Pareto fronts for problems solved using approach III, see Table 6. (results of 10 runs)

V. Conclusions

In this work, we propose a data-driven approach that deals with mixed variables types, continuous and categorical. The proposed approach was applied to an experimental mixed database related to hybrid discontinuous composite materials. Our first goal, in this work, was to build surrogate GP models. The second goal was to make a final recommendation using three different mixed categorical-continuous BO optimization approaches: mono-objective, targeted and multi-objective. The efficiency of our proposed methodology was tested. In particular, the obtained surrogate models were shown to be sufficiently accurate and our optimization process was able to identify correctly optimal fibres.

In this work, we constructed first the surrogate models and then performed the optimization of such models. It is possible to improve the optimization results by feeding directly the database to the optimization approach. Using


Adrien Hakimi
email: adrien.hakimi.etu@univ-lille.fr

Cyrille Bergoin

Patrick Mucci

Adrien Orcid

Hakimi

Evidence of ventilatory constraints during exercise in hypermobile Ehlers-Danlos syndrome

Keywords: breathing, dynamic hyperinflation, expiratory flow limitation, flow-volume loops, hypermobility, lung

Purpose: Hypermobile Ehlers-Danlos syndrome (hEDS) is a connective tissue disorder with many different symptoms such as pain, fatigue, dysautonomia or respiratory symptoms. Among the respiratory manifestations described, the most frequent are exertional dyspnea and breathing difficulties. Mechanical ventilatory constraints during exercise could participate in these respiratory manifestations. The objective of this study was to explore the response of pulmonary flow-volume loops to exercise in patients with hEDS and to look for dynamic hyperinflation and expiratory flow limitation during exercise. Methods: For this purpose, breathing pattern and tidal exercise flow-volume loops were recorded at two workloads (30% and 80% of the peak power output) of a constant load exercise test. Results: Twelve patients were included (11 women, mean age 41±14 years). The results showed a decrease (p = 0.028) in the inspiratory capacity (from 3.12±0.49 L to 2.97±0.52 L), an increase (p = 0.025) in the end-expiratory lung volume (from 0.73±0.68L to 0.88±0.66L, i.e. from EELV comprising 17±12% to 21±12% of forced vital capacity) between the two workloads in favor of dynamic hyperinflation, and half of the patients had expiratory flow limitations. Conclusion: This exploratory study provides evidence for mechanical ventilatory constraints during exercise in patients with hEDS, which may induce discomfort during exercise and could contribute to the respiratory symptomatology.

Main text Introduction

Ehlers-Danlos syndromes (EDSs) are a group of inherited connective tissue disorders mainly characterized by joint hypermobility, skin hyperextensibility, and tissue fragility [START_REF] Malfait | The 2017 international classification of the Ehlers-Danlos syndromes[END_REF]. These syndromes present great clinical and genetic heterogeneity and are classified in 13 subtypes [START_REF] Malfait | The 2017 international classification of the Ehlers-Danlos syndromes[END_REF]. The most common type is hypermobile Ehlers-Danlos syndrome (hEDS), whose diagnosis remains clinical [START_REF] Sulli | Ehlers-Danlos syndromes : State of the art on clinical practice guidelines[END_REF]. In addition to the musculoskeletal manifestations, general symptoms, such as pain and fatigue, are usually presents in patients with hEDS, as well as dysautonomia, gastrointestinal disorders, anxiety, or respiratory symptoms [START_REF] Castori | Natural history and manifestations of the hypermobility type Ehlers-Danlos syndrome : A pilot study on 21 patients[END_REF][START_REF] De Wandele | Clinical heterogeneity in patients with the hypermobility type of Ehlers-Danlos Syndrome[END_REF][START_REF] Tinkle | Hypermobile Ehlers-Danlos syndrome (a.k.a. Ehlers-Danlos syndrome Type III and Ehlers-Danlos syndrome hypermobility type) : Clinical description and natural history[END_REF]. Breathing difficulties and exertional dyspnea are frequently reported [START_REF] Bascom | Respiratory manifestations in the EHLERS-DANLOS syndromes[END_REF][START_REF] Castori | Natural history and manifestations of the hypermobility type Ehlers-Danlos syndrome : A pilot study on 21 patients[END_REF][START_REF] Chohan | A review of respiratory manifestations and their management in Ehlers-Danlos syndromes and hypermobility spectrum disorders[END_REF][START_REF] Hamonet | Manifestations respiratoires dans le syndrome d'Ehlers-Danlos (SED). Nouveaux apports thérapeutiques[END_REF][START_REF] Morgan | Asthma and airways collapse in two heritable disorders of connective tissue[END_REF], but only a few experimental studies have been conducted on lung function in patients with hEDS [START_REF] Reychler | Inspiratory muscle strength training improves lung function in patients with the hypermobile Ehlers-Danlos syndrome : A randomized controlled trial[END_REF].

The respiratory symptoms of hEDS are still not clearly explained and raise many questions [START_REF] Bascom | Respiratory manifestations in the EHLERS-DANLOS syndromes[END_REF].

In patients with hEDS, some studies have reported preserved pulmonary function [START_REF] Reychler | Inspiratory muscle strength training improves lung function in patients with the hypermobile Ehlers-Danlos syndrome : A randomized controlled trial[END_REF], while other studies have found increased lung volumes [START_REF] Morgan | Asthma and airways collapse in two heritable disorders of connective tissue[END_REF], restrictive or mixed respiratory insufficiencies [START_REF] Castori | Natural history and manifestations of the hypermobility type Ehlers-Danlos syndrome : A pilot study on 21 patients[END_REF], gas transfer anomalies [START_REF] Ayres | Abnormalities of the lungs and thoracic cage in the Ehlers-Danlos syndrome[END_REF][START_REF] Morgan | Asthma and airways collapse in two heritable disorders of connective tissue[END_REF] or deficits in inspiratory muscles strength [START_REF] Reychler | Inspiratory muscle strength training improves lung function in patients with the hypermobile Ehlers-Danlos syndrome : A randomized controlled trial[END_REF]. Anomalies in lung or chest wall compliance and deformities of the thoracic cage have also been reported [START_REF] Ayres | Abnormalities of the lungs and thoracic cage in the Ehlers-Danlos syndrome[END_REF][START_REF] Morgan | Asthma and airways collapse in two heritable disorders of connective tissue[END_REF]. These last impairments may be at the origin of mechanical ventilatory constraints during exercise, leading to dyspnea [START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF].

During exercise, in order to respond to the increase in metabolic demand, increased ventilation leads to an expansion of the tidal volume. In healthy young subjects, this tidal volume response results in the utilization of inspiratory and expiratory reserves (i.e., to a decrease in end-expiratory lung volume [EELV] and an increase in end-inspiratory lung volume [EILV]) [START_REF] Milne | Evaluation of Dynamic Respiratory Mechanical Abnormalities During Conventional CPET[END_REF]. However, in the case of an increase in lung compliance (e.g., due to collagen or extracellular matrix abnormalities that can be found in EDS), the energy stored at the end of inspiration could be reduced and lead to an early involvement of the expiratory muscles and difficulties in ensuring sufficient exhalation. In this case, the increase in tidal volume is achieved with a temporary increase in EELV above its resting level and a decrease in inspiratory capacity [START_REF] Langer | Lung hyperinflation in chronic obstructive pulmonary disease : Mechanisms, clinical implications and treatment[END_REF]. This breathing pattern is called dynamic hyperinflation and leads to early ventilatory limitation and functional weakness of the inspiratory muscles [START_REF] Langer | Lung hyperinflation in chronic obstructive pulmonary disease : Mechanisms, clinical implications and treatment[END_REF]. It may be hypothesized that mechanical ventilatory constraints, such as dynamic hyperinflation and/or expiratory flow limitation, could exist in patients with hEDS during exercise.

To our knowledge, the pulmonary mechanical constraints during exercise have not been studied in patients with hEDS. These constraints are commonly explored using breathing reserve, which is a comparison of the ventilation achieved during exercise and the maximal voluntary ventilation (MVV) measured at rest [START_REF] Guenette | Inspiratory Capacity during Exercise : Measurement, Analysis, and Interpretation[END_REF][START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF][START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF][START_REF] Nourry | Evidence of Ventilatory Constraints in Healthy Exercising Prepubescent Children[END_REF]. But, breathing reserve does not provide precise information on the type or source of ventilatory constraints [START_REF] Guenette | Inspiratory Capacity during Exercise : Measurement, Analysis, and Interpretation[END_REF][START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF][START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF][START_REF] Nourry | Evidence of Ventilatory Constraints in Healthy Exercising Prepubescent Children[END_REF]. It was shown that the use of tidal exercise flow-volume loops plotted within the maximal flow-volume loop provide more specific information, in particular about the breathing strategy during exercise [START_REF] Guenette | Inspiratory Capacity during Exercise : Measurement, Analysis, and Interpretation[END_REF][START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF][START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF][START_REF] Nourry | Evidence of Ventilatory Constraints in Healthy Exercising Prepubescent Children[END_REF]. Furthermore, assessment of the response of flow-volume loops to exercise is a useful tool for the investigation of exercise intolerance and exertional dyspnea [START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF][START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF]. The degree of expiratory flow limitation (EFL) is an index of the balance between ventilatory demand and ventilatory capacity combined with how the subject chooses to regulate their EELV [START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF].

The objective of this study was to explore the mechanical ventilatory constraints to exercise and to investigate the response of flow-volume loops at two exercise intensities in order to determine the achievement of dynamic hyperinflation and/or EFL in patients with hEDS. The hypothesis was that patients with hEDS could present ventilatory constraints to exercise as expiratory flow limitations or hyperinflation.

Materials and Methods

Patients

The participants were recruited among the patients who were undergoing a rehabilitation program in the Clinique de la Mitterie (Lomme, France). The rehabilitation program is based on physical activity (2/3) and on educational or mental well-being activities (1/3). It was designed for hEDS patients and is described elsewhere [START_REF] Hakimi | Immediate and 6-week after effects of a rehabilitation program for EHLERS-DANLOS syndrome hypermobile type patients : A retrospective study[END_REF]. Inclusion criteria were having performed recent cardiopulmonary exercise testing and having been diagnosed with hEDS according to the 2017 criteria [START_REF] Malfait | The 2017 international classification of the Ehlers-Danlos syndromes[END_REF]. All patients were volunteers and gave their informed consent. This study is part of a clinical trial that has been approved by a national ethics committee (CPP Ile de France IV, 2019/71).

Study design

First, patients performed an incremental cardiopulmonary exercise test before their inclusion in a rehabilitation program. Then, during the first week of the rehabilitation program, patients performed spirometry, followed by a constant load exercise test (CLET) with measurement of pulmonary flow-volume loops.

Incremental cardiopulmonary exercise test

This incremental test was conducted on a calibrated, electromagnetically braked cycle ergometer (Ergoselect 200, Ergoline, Germany), with a computerized open-circuit gas collection system (Ergocard, Medisoft, Belgium).

Before the start of each test, the airflow sensor was calibrated with a 3 L syringe, and gas analyzers were calibrated against known gases. It consisted of a 1-minute warm-up (the warm-up workload varied according to the fitness level), followed by a progressive increase of workload in a ramp of 10 W/min. The peak power output (PPO) was determined during this test. Heart rate (HR) and cycle by cycle data, including oxygen consumption (V ̇O2), carbon dioxide output (V ̇CO2), breathing frequency (fb), tidal volume (VT), and minute ventilation (V ̇E), were measured.

Predicted maximal HR (HRmax, pred) was estimated as follows: HRmax, pred = (208-0.7) x age [START_REF] Tanaka | Age-predicted maximal heart rate revisited[END_REF].

The respiratory exchange ratio (RER) was calculated as follow: RER = V ̇CO2/ V ̇O2. Predicted maximal oxygen consumption was estimated as recommended by the American Thoracic Society with Wasserman equations (American Thoracic Society & American College of Chest Physicians, 2003;Wasserman, 2005).

Spirometry

A forced spirometry maneuver was conducted just before the CLET following the current recommendations [START_REF] Graham | Standardization of Spirometry 2019 Update. An Official American Thoracic Society and European Respiratory Society Technical Statement[END_REF]. The material (the same as for the incremental test) was calibrated before each test. The spirometry allowed us to establish a maximal pulmonary flow-volume envelope. No bronchodilators were used.

The equations used to calculate the predicted values are the global lung function initiative 2012 equations [START_REF] Quanjer | Multi-ethnic reference values for spirometry for the 3-95-yr age range : The global lung function 2012 equations[END_REF]. The MVV was estimated as follows: MVV = FEV1 x 40 [START_REF] Campbell | A comparison of the maximum voluntary ventilation with the forced expiratory volume in one second : An assessment of subject cooperation[END_REF], where FEV1 is the forced expiratory volume in 1 second measured with spirometry.

Constant load exercise test (CLET)

The CLET was conducted on the same ergometer with the same measurement system as the incremental cardiopulmonary exercise test. Patients started with 3 minutes of pedaling at a moderate workload of 30% of PPO. This was followed by 2 minutes of linear increase to 80% of PPO. Then, the patient was instructed to pedal for as long as possible at a target workload corresponding to 80% of PPO, with a maximum of 20 minutes. The target pedaling frequency was 60 rpm. If the patient could no longer maintain the pedaling frequency (55 rpm minimum) or reach the 20 minute limit, the test was stopped and the patient began a recovery period at 30% of PPO. The same data as for the incremental test were recorded, in addition to inspiratory time (Ti), total respiratory time (Ttot) and arterial oxygen saturation as indicated by pulse oximetry (SpO2). Dyspnea was assessed at rest before and at the end of the exercise testing by the modified Borg scale [START_REF] Borg | Psychophysical bases of perceived exertion[END_REF].

In order to investigate ventilatory constraints, we performed tidal exercise flow-volume loops by using inspiratory capacity (IC) maneuvers at the third minute of the exercise at 30% of PPO and at the end of the exercise at 80% of PPO. These maneuvers correct the relative position of the tidal volume in the maximal flow-volume envelope during exercise and permit us to calculate the EELV and EILV in percentages of the forced vital capacity (FVC) as follows [START_REF] Guenette | Inspiratory Capacity during Exercise : Measurement, Analysis, and Interpretation[END_REF]:

EELV = [(FVC-IC) x 100] / FVC EILV = [(EELV + VT) x 100] / FVC.
Dynamic hyperinflation is considered in this study when EELV at 80% of PPO equals or exceeds the 30% of PPO value [START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF]. The EFL was graphically calculated as the volume that meets or exceeds the expiratory boundary of the maximal flow-volume loop and was expressed as a percent of VT as previously described [START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF].

Breathing reserve was estimated from MVV and V ̇E reached during exercise, according to the following equation [START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF]:

BR (%) = [(MVV-V ̇E) / MVV] x 100.
Ventilatory efficiency, which represents the quantity of ventilation needed to eliminate metabolically produced CO2, was evaluated using the V ̇E/V ̇CO2 ratio [START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF][START_REF] Sun | Ventilatory Efficiency during Exercise in Healthy Subjects[END_REF]. A mismatching between ventilation and CO2 production can contribute to exercise dyspnea [START_REF] Sun | Ventilatory Efficiency during Exercise in Healthy Subjects[END_REF].

All cycle by cycle data measured during the CLET were averaged by minute. Every minute where an IC maneuver was conducted was removed from the analysis. The 3 minutes of exercise at 30% of PPO were averaged and used as the 30% of PPO values. The last minute of the exercise at 80% of PPO without IC maneuvers was used as the 80% of PPO value. If a subject failed to achieve at least 1 minute at 80% of PPO, the averaged minute was replaced by the last value measured at 80% of PPO.

Statistical analysis

Data are described with mean and standard deviation (SD). Statistical analysis was performed on SigmaStat Version 3.5 (Systat Software Inc). Before each test, normality of the data distribution and equality of variances were tested. Differences between the 30% of PPO values and the 80% of PPO values from the CLET were tested with a paired t-test or a Wilcoxon signed rank test. Results were considered as significant for a p-value < 0.05.

Results

Patients

Twelve patients, including 11 women and 1 man, were included in this study. All patients were diagnosed with hEDS according to the 2017 criteria. Anthropometric data are presented in Table 1.

Incremental cardiopulmonary exercise test

Data from the incremental test are presented in Table 2. Overall, the patients showed limited aerobic fitness, with a mean of 79% of the predicted V ̇O2peak. Only two patients conserved a V ̇O2peak of more than 100% of the predicted value. All other patients were between 57% and 83% of the predicted V ̇O2peak values.

Spirometry

The spirometry data were within the normal range (Table 1), except for one patient with FEV1 at 74% of the predicted value and FVC at 73% of the predicted value but a FEV1/FVC ratio of 83% (102% of the predicted value). Three subjects present FEV1 values ≥ 100% of the predicted value, four subjects between 90% and 99% and four subjects between 80% and 89%. In the same way, two subjects had FVC values ≥ 100% of the predicted value, six subjects between 90% and 99% and three subject between 80% and 89%. The FEV1/FVC ratios were between 75% and 93%.

Constant load exercise test

No patient reached the maximum 20 minutes on the cycle exercise at 80% of PPO. The average time to exhaustion at 80% of PPO was 200 ± 165 seconds. Most of the patients (n=8) stopped for exhaustion with inability to further maintain the pedaling frequency of 55 rpm, 2 patients stopped for dyspnea, and 2 patients because of articular pain in the legs. The mean cardiorespiratory data and the breathing pattern corresponding to 30% of PPO and 80% of PPO are presented respectively in Tables 3 and4. There was a significant (p < 0.05) increase in the V ̇O2, HR, V ̇E, fR, VT and Ti/Ttot ratio. There was a significant increase in dyspnea assessed by the modified Borg scale from 0.9 ± 1.4 at rest to 5.5 ± 1.7 at the end of the test. No significant change for neither the SpO2 nor for the V ̇E/V ̇CO2 ratio was observed.

Regarding the mechanical ventilatory constraints to exercise, breathing reserve significantly (p < 0.001) decreased from 83% ± 4% to 60% ± 12% between 30% and 80% of PPO. A significant decrease in the IC (p = 0.028) and a significant increase in EELV (p = 0.025) associated with a significant increase in EILV (p < 0.001) were observed (Table 5 and Figure 1). All patients showed an increase in EILV (from 4.0% to 25.3% of FVC). Three patients showed a decrease in EELV (4.6%, 3.8% and 0.2% of FVC, respectively) and nine patients showed dynamic hyperinflation (i.e., an increase from 1.3% to 13.4% of FVC). In addition, EFLs were found in six patients (50%), with a mean limitation equal to 46% of VT (from 16% to 75%) at 30% of PPO and equal to 43% of VT (from 33% to 58%) at 80% of PPO. There was no significant change for the EFL between the exercise at 30% of PPO and the exercise at 80% of PPO.

Discussion

The main objective of this study was to explore the pulmonary mechanical constraints to exercise in patients with hEDS and to investigate the response of flow-volume loops to exercise. The results showed a significant increase in EELV during exercise in favor of dynamic hyperinflation, and half of the patients had EFL.

The majority of patients showed reduced peak oxygen consumption in comparison with predicted values during the incremental test. This is an indicator of reduced exercise capacity (American Thoracic Society & American College of Chest Physicians, 2003). As mentioned by some authors, fear of trauma and pain can lead to deconditioning and exercise intolerance, which could contribute to the exertional dyspnea observed during the CLET [START_REF] Baeza-Velasco | Cognitive, emotional, and behavioral considerations for chronic pain management in the Ehlers-Danlos syndrome hypermobilitytype : A narrative review[END_REF][START_REF] Syx | Hypermobility, the Ehlers-Danlos syndromes and chronic pain[END_REF]. The cycling time to exhaustion recorded by the patients at 80% of PPO (200 ± 165 seconds) was relatively low. Indeed, it is comparable to the times observed in patients with chronic obstructive pulmonary disease (about 4 minutes) at 75% of their PPO (van 't Hul et al., 2003). This supports impaired exercise capacity in patients with hEDS. However, we have to take into account that the participants in this study were recruited among future participants in a rehabilitation program and that respiratory comorbidities were not recorded although this could have been interesting.

At the end of the incremental exercise or during CLET, breathing reserves were relatively substantial (around 60% during CLET at 80% PPO), as it is generally reported in healthy subjects (American Thoracic Society & American College of Chest Physicians, 2003). This result did not support existence of mechanical ventilatory constraints during exercise. Nevertheless, it was not surprising because, during CLET, the exercise intensities were submaximal and the spirometric data were mostly within the normal ranges. However, analysis of breathing reserve alone is insufficient to estimate mechanical limitations of breathing, as largely described in the literature [START_REF] Guenette | Inspiratory Capacity during Exercise : Measurement, Analysis, and Interpretation[END_REF][START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF][START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF][START_REF] Nourry | Evidence of Ventilatory Constraints in Healthy Exercising Prepubescent Children[END_REF]. The use of tidal exercise flow-volume loops during exercise has provided additional information on the respiratory pattern of these patients. The significant increase in EELV found during exercise is an argument for dynamic hyperinflation in patients with hEDS. For these patients, the tidal volume expansion during exercise is achieved by a major reduction in the inspiratory reserve volume (represented by the increase in EILV), which is normally compensated in healthy subjects by a decrease in EELV [START_REF] Milne | Evaluation of Dynamic Respiratory Mechanical Abnormalities During Conventional CPET[END_REF]. This theoretical reduction in EELV allows the optimal positioning of the VT on the pressure-volume curve of the respiratory system [START_REF] Langer | Lung hyperinflation in chronic obstructive pulmonary disease : Mechanisms, clinical implications and treatment[END_REF][START_REF] Milne | Evaluation of Dynamic Respiratory Mechanical Abnormalities During Conventional CPET[END_REF]. As there is no decrease here, the reduction in length of the inspiratory muscles translates to a diminished force generating capacity and to an increased elastic load [START_REF] Langer | Lung hyperinflation in chronic obstructive pulmonary disease : Mechanisms, clinical implications and treatment[END_REF].

In addition, EFL was observed in half of the patients. It may be due to abnormalities such as the tendency toward airway collapse, as highlighted by [START_REF] Morgan | Asthma and airways collapse in two heritable disorders of connective tissue[END_REF]. This tendency toward EFL may be one of the causes of dynamic hyperinflation, with the need to increase operational lung volumes to reduce these flow limitations.

However, as half of the patients do not have EFL, this hypothesis does not fully explain the general trend towards hyperinflation.

Based on the hypothesis that lung compliance may be increased in hEDS [START_REF] Morgan | Asthma and airways collapse in two heritable disorders of connective tissue[END_REF], the energy stored at the end of inspiration could be reduced. In this case, active and premature involvement of the expiratory muscles may be necessary. This early work of the expiratory muscles associated with a potential limited strength of these muscles [START_REF] Reychler | Inspiratory muscle strength training improves lung function in patients with the hypermobile Ehlers-Danlos syndrome : A randomized controlled trial[END_REF] could explain the absence of a decrease in EELV. However, this mechanism should be compensated by the drop in compliance as EILV approaches the maximal lung volume [START_REF] Johnson | Emerging Concepts in the Evaluation of Ventilatory Limitation During Exercise[END_REF].

Therefore, compliance impairment does not completely explain hyperinflation either. Further investigations are needed in order to clarify the involved mechanisms. A study with electromyographical activation of respiratory muscles during exercise and of lung compliance in hEDS could provide further information.

In addition, the inspiratory muscle weakness that has been shown in patients with hEDS [START_REF] Reychler | Inspiratory muscle strength training improves lung function in patients with the hypermobile Ehlers-Danlos syndrome : A randomized controlled trial[END_REF] could be also at the origin of a decrease in the velocity of inspiratory muscle contraction, thus an increase in inspiratory time leading to dynamic hyperinflation. In our case, the inspiratory time increase between 30% and 80% of PPO from 39% of the respiratory total time to 45%. Although these values increase, they remain below or within the described values from healthy sedentary individuals of between 45% and 47% [START_REF] Neder | The pattern and timing of breathing during incremental exercise : A normative study[END_REF]. The changes in inspiratory times alone, therefore, do not explain the dynamic hyperinflation.

Proprioceptive impairment has been shown in hEDS [START_REF] Rombaut | Joint position sense and vibratory perception sense in patients with Ehlers-Danlos syndrome type III (hypermobility type)[END_REF][START_REF] Scheper | The association between muscle strength and activity limitations in patients with the hypermobility type of Ehlers-Danlos syndrome : The impact of proprioception[END_REF]. It is possible that these sensory disorders also affect the proprioceptive receptors of the rib cage and respiratory muscles. Thus, poor perception of lung volumes could also be a partial explanation for the dynamic hyperinflation. Abnormalities in sensory perception could also influence the thresholds of dyspnea perception, particularly in a predominantly female population in which greater dyspnea perception is found compared to men [START_REF] Archiza | An integrative approach to the pulmonary physiology of exercise : When does biological sex matter[END_REF]. However, these hypothesis need to be explored in further studies. Finally, as frequently described in the literature [START_REF] Chohan | A review of respiratory manifestations and their management in Ehlers-Danlos syndromes and hypermobility spectrum disorders[END_REF], the results showed significant dyspnea in patients with hEDS during the CLET. The observed V ̇E/V ̇CO2 ratios of 38.4% at 30% of PPO and of 37.5% at 80% of PPO are above the recommended value of 34% for establishing ventilatory inefficiency [START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF]. This ventilatory inefficiency associated with the dynamic hyperinflation may reflect neuromotor uncoupling and could partly explain the exertional dyspnea [START_REF] Neder | The role of evaluating inspiratory constraints and ventilatory inefficiency in the investigation of dyspnea of unclear etiology[END_REF].

The exploratory nature of this study is associated with certain methodological limitations such as the absence of measurement of inspiratory capacity at rest, which would have provided information on the presence of possible hyperinflation at rest, and the absence of measurement of the Borg scale at 30% PPO. Future studies should investigate mechanical ventilatory constraints in hEDS more systematically with repeated measurements of inspiratory capacity, including at rest, as well as dyspnea.

Conclusion

In conclusion, this exploratory study provides, for the first time, evidence for ventilatory constraints during exercise in patients with hEDS, such as dynamic hyperinflation and EFL. These constraints may induce discomfort during exercise. Further investigations of lung function in hEDS, during exercise and at rest, are needed to better understand the underlying mechanisms and to generalize to a wider patient population with hEDS. van 't Hul, A., Gosselink, R., & Kwakkel, G. (2003)

Statements and declarations

FiguresFig. 1

 1 Figures

  Fig. 1 Individual changes in EELV and EILV for all patients (n = 12) between 30% of PPO and 80% of PPO.The male participant is represented by empty markers.

  . Constant-load cycle endurance performance : Test-retest reliability and validity in patients with COPD. Journal of Cardiopulmonary Rehabilitation, 23(2),

	143-150. https://doi.org/10.1097/00008483-200303000-00012
	Wasserman, K. (Éd.). (2005). Principles of exercise testing and interpretation : Including pathophysiology and
	clinical applications (4th ed). Lippincott Williams & Wilkins.

Table 1 :

 1 Anthropometric and spirometry data.

		Mean ± SD
	Data (n = 12)	[Min : Max]
		41 ± 14
	Age (years)
		[23 : 65]
		172 ± 9
	Height (cm)	[159 : 192]
		79 ± 13
	Weight (kg)
		[58 : 99]
	BMI (kg/m 2 )	26.6 ± 3.6
		[21.4 : 32.0]
		5.2 ± 1.5
	Beighton score
		[3 : 7]
	FVC (L)	3.85 ± 0.84
	(% predicted value)	(92 ± 9%)
	FEV1 (L)	3.18 ± 0.75
	(% predicted value)	(93 ± 11%)
	FEV1/FVC (%)	82 ± 5
	(% predicted value)	(100 ± 6%)
	MVV (L/min)	127 ± 30

Abbreviations: BMI, body mass index; FVC, forced vital capacity; FEV1, forced expiratory volume in 1 second; MVV, maximal voluntary ventilation.

Table 2 :

 2 Ventilatory threshold and peak data during incremental cardiopulmonary exercise.

	Data (n = 12)

Abbreviations: V ̇O2peak, oxygen uptake at peak exercise; HR, heart rate; V ̇E, minute ventilation; MVV, maximal voluntary ventilation; fb, breathing frequency; VT, tidal volume.

Table 3 :

 3 Mean cardiorespiratory data during constant load exercise test.Values are mean ± SD. Abbreviations: V ̇O2, oxygen uptake; HR, heart rate; SpO2, arterial oxygen saturation.

	Data (n = 12)	30% PPO	80% PPO
	Workload (W)	32 ± 10	86 ± 27*
	V ̇O2 (L/min)	0.67 ± 0.12	1.33 ± 0.33*
	(% predicted value)	(38 ± 10)	(76 ± 25)
	HR (bpm)	102 ± 14	142 ± 23*
	(% predicted value)	(57 ± 8)	(79 ± 14)
	SpO2 (%)	96.4 ± 1.8	97.1 ± 1.5

*Significant difference with 30% PPO.

Table 4 :

 4 Breathing pattern during constant load exercise test.Values are mean ± SD. Abbreviations: V ̇E, minute ventilation; MVV, maximal voluntary ventilation; fb, breathing frequency; VT, tidal volume; V ̇CO2, carbon dioxide output; Ti, inspiratory time; Ttot, total respiratory time.

	Data (n = 12)	30% PPO	80% PPO
	V ̇E (L/min)	20.5 ± 4.4	49.8 ± 16.4*
	Breathing reserve (%)	83 ± 4	60 ± 12*
	fb (breaths/min)	20 ± 4	26 ± 6*
	VT (L)	1.08 ± 0.22	1.91 ± 0.30*
	V ̇E/V ̇CO2	38.4 ± 3.9	37.5 ± 7.9
	Ti/Ttot	0.39 ± 0.04	0.45 ± 0.04*

*Significant difference with 30% PPO.

Table 5 :

 5 Pulmonary flow volume data during constant load exercise.Values are mean ± SD, except for number of subjects with EFL which is the number of subjects and the percentage of the total number of subjects.

	30% PPO	80% PPO

Abbreviations: IC, inspiratory capacity; EELV, end-expiratory lung volume; FVC, forced vital capacity; EILV, end-inspiratory lung volume; VT, tidal volume; EFL, expiratory flow limitation. IC, EELV, and EILV are expressed for all the subjects (n = 12). EFL in L and %VT are expressed only for the six expiratory flow limited subjects. *Significant difference with 30% PPO.

Funding

The authors did not receive support from any organization for the submitted work.

Data availability statement

The data that support the findings of this study are available from the corresponding author upon reasonable request.

Competing interests

The authors have no relevant financial or non-financial interests to disclose.

Compliance with Ethical Standards:

This study is part of a clinical trial that has been approved by a national ethics committee (CPP Ile de France IV, 2019/71). All patients were volunteers and gave their informed consent.

Authors' contribution

AH, CB, and PM contributed to the design and implementation of the research, to the analysis of the results and to the writing of the manuscript.