Siyang Deng

Reda El Bechari

Stéphane Brisset

Stéphane Clénet

Iterative Kriging-based Methods for Expensive Black-Box Models

Reliability-Based Design Optimization (RBDO) in electromagnetic field problems requires the calculation of probability of failure leading to a huge computational cost in the case of expensive models. Three different RBDO approaches using kriging surrogate model are proposed to overcome this difficulty by introducing an approximation of the objective function and constraints. These methods use different infill sampling criteria (ISC) to add samples in the process of optimization or/and in the reliability analysis. Several enrichment criteria and strategies are compared in terms of number of evaluations and accuracy of the solution.

Index Terms-Infill sampling criteria, kriging model, reliability analysis, reliability-based design optimization.

I. INTRODUCTION ELIABILITY-BASED DESIGN OPTIMIZATION

(RBDO) approaches can be divided into Double-Loop (DLM), Single-Loop (SLM) and Sequential Decoupled Methods (SDM). They have emerged in the past few decades and become more and more popular in electromagnetics owing to their ability to account for uncertain parameters. However, for expensive black-box models, the computational burden can become unbearable.

To overcome this issue, iterative kriging surrogate models have been proposed to reduce the number of evaluations [START_REF] Lee | A sampling technique enhancing accuracy and efficiency of metamodel-based RBDO: Constraint boundary sampling[END_REF]- [START_REF] Lee | Sampling-based RBDO using the stochastic sensitivity analysis and Dynamic Kriging method[END_REF]. Infill Sampling Criterion (ISC) was used with the aim of improving the quality of the surrogate model, and searching for the solution of the optimization problem. However, the metamodel is established before starting RBDO and no enrichment is made during neither optimization nor reliability analysis.

With the purpose of enhancing the efficiency, different strategies including the choice of the ISC and the positioning of sample enrichments in the optimization process are investigated in this paper for each aforementioned type of RBDO approaches, so that the reliabilities are also analyzed by metamodels. A mathematical example is used to compare with classic RBDO, i.e. without kriging model, and highlight the most effective strategy. Then, RBDO of a transformer modelled by a time consuming model based on the Finite Element method is performed with the most effective strategy.

II. INFILL SAMPLING CRITERIA

Iterative surrogate-based optimization methods start with a small set of initial sampling points to create a preliminary metamodel. Then, the infill sampling criteria are considered as new objective functions to add points into the sample set and update the meta-model until the predicted error is less than a chosen tolerance. A great advantage of this approach is that it enhances the accuracy of meta-model and search the probabilistic optimum simultaneously with a small amount of samples.

Expected Improvement (EI) criterion [START_REF] Jones | Efficient global optimization of expensive black-box functions[END_REF] is widely used for surrogate-based optimizations without constraints.

𝐸𝐼 𝑓 = { (𝑓 𝑚𝑖𝑛 -𝑓 ̂)𝛷(𝑧) + 𝑠̂𝑓𝜙(𝑧) if 𝑠̂𝑓 > 0 0 if 𝑠̂𝑓 = 0 (1)
where 𝑓 𝑚𝑖𝑛 is the best current sampled objective function value, 𝑓 ̂ and 𝑠̂𝑓 are the predicted value and the mean square error (MSE), 𝜙(•) and 𝛷(•) denote the probability density function and the cumulative distribution function of the standard normal distribution respectively, and 𝑧 = (𝑓 𝑚𝑖𝑛 -𝑓 ̂) 𝑠̂𝑓 ⁄ . However, as EI is multimodal, more attention should be paid on the infill criterion to be sure to find the global solution. The Weighted EI (WEI) criterion [START_REF] Xiao | Adaptive weighted expected improvement with rewards approach in Kriging assisted electromagnetic design[END_REF] seems to be more suitable as it adds weights into EI expression to balance exploration (right part) and intensification (left part).

W𝐸𝐼 𝑓 = { 𝜔(𝑓 𝑚𝑖𝑛 -𝑓 ̂)𝛷(𝑧) + (1 -𝜔)𝑠̂𝑓𝜙(𝑧) if 𝑠̂𝑓 > 0 0 if 𝑠̂𝑓 = 0 (2)
Choosing small weight 𝜔 prevents WEI from converging to a local minimum if the initial sampling is inside the security domain. This condition is quite difficult to satisfy for many devices, as their security domains may be small and sometimes discontinuous. To avoid this issue, a Modified WEI (MWEI) combined with the surrogate objective function is proposed:

MW𝐸𝐼 𝑓 = W𝐸𝐼 𝑓 -𝜔𝑓 ̂ (3)
Investigations on the same example as in [START_REF] Xiao | Adaptive weighted expected improvement with rewards approach in Kriging assisted electromagnetic design[END_REF] show that a weight equal to 0.1 provide to a global optimum with less iterations.

For constrained problems, an extended method consists in multiplying the value of EI by the probability of feasibility (PF) [START_REF] Schonlau | Computer Experiments and Global Optimization[END_REF]. However, PF may prevent the sampling on the constraint boundary where the deterministic optimum may lie. Another constraint handling method is the Expected Violation (EV) method [START_REF] Audet | A surrogate-model-based method for constrained optimization[END_REF] but the number of candidate points to evaluate can be very large. An alternative method is to use the predicted value of the constraint functions 𝑔 ̂ directly as constraints in the infill sub-problem [START_REF] Sasena | Flexibility and efficiency enhancements for constrained global design optimization with kriging approximations[END_REF].

R

Manuscript received June, 28 2017; revised September 1, 2017; accepted. Date of publication; date of current version. Corresponding author: Siyang Deng (e-mail: siyang.deng@centralelille.fr).

III. INFILL STRATEGIES FOR RBDO METHODS

RBDO is a combination of deterministic constrained optimizations and reliability analysis. For the first one, the design variables 𝑑 are the mean value of random variables 𝑋 and the standard deviations 𝜎 are constant. For reliability analysis, 𝑑 and 𝜎 are constants and the design variable 𝑥 is a realization of 𝑋. 𝛽 𝑡 is the given target reliability index.

A. Double-Loop Method

DLM like Performance Measure Approach (PMA) [START_REF] Tu | A new study on reliability-based design optimization[END_REF] has a nested structure: The outer loop seeks for the optimum and the inner loop searches the Most Performance Target Point (MPTP) that maximize the constraint subject to a given reliability index.

There are two places where ISC can be introduced to improve the accuracy of kriging model: Outer loop and inner loop. For the inner loop, EI of 𝑔 is used to find MPTP by solving the optimization problem in Eq. ( 4)-( 5) as the constraint on the reliability index is an explicit function of the design variables:

𝑥 * = argmax 𝑥 𝐸𝐼 𝑔 (𝑥) 𝑠. 𝑡. ‖(𝑥 -𝑑) 𝜎 ⁄ ‖ = 𝛽 𝑡 ( 4
)
𝐸𝐼 𝑔 = { (𝑔 ̂-𝑔 𝑚𝑎𝑥 )𝛷(𝑧 𝑔 ) + 𝑠̂𝑔𝜙(𝑧 𝑔 ) if 𝑠̂𝑔 > 0 0 if 𝑠̂𝑔 = 0 (5)
where 𝑧 𝑔 = (𝑔 ̂-𝑔 𝑚𝑎𝑥 ) 𝑠̂𝑔 ⁄ , 𝑠̂𝑔 is the MSE of the constraint, and 𝑔 𝑚𝑎𝑥 is the maximum sampled constraint value.

For outer loop, the criterion MWEI is preferred to avoid local solutions and the implicit inequality constraints are computed by the inner loop. However, as the two loops are nested, the enrichment in inner loop may bring out thousands of model evaluations. To test it, two strategies are proposed: the first one (PMA1) adds new samples only inside the outer loop, whereas the second (PMA2) enriches inside both loops.

B. Single-Loop Method

For SLM like Single Loop Approach (SLA) [START_REF] Liang | A single-loop method for reliability-based design optimization[END_REF], the main point is that the inner loop optimization is replaced by an approximation based on a first order Taylor expansion to avoid the numerous evaluations required to find the MPTP.

It is important to note that due to its approximation, the method itself has already loose some precision. Therefore, it is expected that with a surrogate model, the two approximation errors will be added and the accuracy will be further reduced.

C. Sequential Decoupled Method

SDM like Sequential Optimization and Reliability Assessment (SORA) [START_REF] Du | Sequential optimization and reliability assessment method for efficient probabilistic design[END_REF] are based on a series of sequential deterministic optimizations and reliability assessments. The main point is to shift the boundaries of constraints inside the feasible domain based on the reliability information obtained in the former iteration. The first optimization aims at searching the global deterministic optimum. Reliability assessment is then conducted to locate the MPTP corresponding to the target reliability index. Finally, new optimizations are carried out by taking into account the shift 𝑡 computed with MPTP.

Three strategies are proposed. In the first one (SORA1), the reliability analysis is the same as in the inner loop of PMA and enrichment is made with EI criterion. For the deterministic optimizations, the constraints are computed with the metamodel 𝑔 ̂(𝑑 -𝑡) and the MWEI is preferred in order to find the global solution.

The second strategy (SORA2) differs from the first one by the fact that enrichment of the kriging models with MWEI criterion is made at first iteration only. For all other iterations, the deterministic optimization is made with the meta-models 𝑓 ̂ and 𝑔 ̂.

For the third strategy (SORA3), if the deterministic optimum found in 𝑘-th cycle 𝑑 𝑘 is close to any optimum of other 𝑘 -1 cycles, as the former reliability assessments have already added points in this region, the accuracy is considered to meet the requirement so there is no need to add samples any more. The proximity criterion defined in ( 6) is checked before entering reliability analysis:

‖(𝑑 𝑘 -𝑑 𝑖 )/𝜎‖ < 𝛽 𝑡 , 𝑖 = 1, … , 𝑘 -1 (6)
where 𝑑 𝑖 is the deterministic optimum found by the 𝑖-th cycle. If ( 6) is satisfied, the meta-model of constraints is used directly and only MPTP are evaluated. For the deterministic optimizations, it takes the same strategy as SORA2. The flowchart of SORA3 is shown in Fig. 1 where 𝑥 𝑘 is the MPTP of 𝑘-th cycle.

IV. COMPARISON OF STRATEGIES

To assess the efficiency of kriging-based RBDO methods, the mathematical example in [START_REF] Kim | A Single-Loop Strategy for Efficient Reliability-Based Electromagnetic Design Optimization[END_REF] with two variables and three constraints is analyzed. Noting that the random variables are Gaussian and their standard deviations are all equal to 0.3. Lower and upper bounds are 0 and 10 respectively for both variables. The target reliability index 𝛽 𝑡 is chosen equal to 2, so that the target probability of failure 𝑃 𝑡 = 𝛷(-𝛽 𝑡 ) = 2.28%. The results are given in Table 1 with an initial sampling of 20 points. The probability of failure 𝑃 𝑓 is calculated by Monte-Carlo Simulation (MCS) with 10 6 samples. For comparison purpose, results given by classic RBDO methods without kriging are also presented.

All the iterative kriging-based RBDO methods lead to a reduced number of evaluations. SLA with kriging has the minimum number of evaluations but is not accurate enough as the maximum probability of failure is much greater than 𝑃 𝑡 due to the approximation used to simplify the reliability analysis. As expected, PMA with infill during inner loops requires thousands of samples to evaluate. The other PMA strategy is faster but the accuracy is not sufficient. Kriging-based SORA strategies lead to the best result and the third one SORA3 is the most efficient. Fig. 2 shows two iterations of SORA3.

As SORA3 seems to be the most efficient meta-model strategy on this mathematical example, it is tested on the RBDO of a transformer with FEM.

V. ELECTROMAGNETIC DEVICE

The electromagnetic device is a single-phase safety isolating transformer with grain-oriented E-I laminations designed for installation in electric cabinet [START_REF] Tran | A Benchmark for Multi-Objective, Multi-Level and Combinatorial Optimizations of a Safety Isolating Transformer[END_REF]. The primary and secondary windings are wound around the frame surrounding the central core (Fig. 3, left).

A. Finite Element Models

Thermal and magnetic phenomena are modeled by using 3D FEA on the eighth of transformer due to symmetries. There are about 43,000 nodes and 290,000 edges in the model. The right part of Fig. 3 shows the mesh in the magnetic circuit, the insulating, the air gap, the frame and the opposing direction of currents in the primary and secondary windings that create flux in the gap between the coils (leakage flux).

For the electromagnetic modeling, all magnetic and electric quantities are assumed sinusoidal. Full-load and no-load simulations are used to compute all the characteristics. The iron losses are computed with the Steinmetz formula and the leakage inductances are calculated with the magnetic co-energy. The core magnetic nonlinearity is taken into account.

In the thermal modeling, some assumptions are considered: the insulator between the core and the coils is in perfect contact with both parts; there is no thermal contact between the exterior coil and the magnetic circuit; there is no thermal exchange with the air trapped between the coils and the iron; there is no convection on the upper and lower sides of the coil; there is no temperature gradient in the copper and the iron, and all surfaces have the same convection coefficient.

A magneto-thermal weak coupling is considered and the computational time is equal to10 minutes on a single core of an Intel Xeon CPU E5-2690 at 2.60 GHz. The copper and iron losses are computed with the magnetic AC solver and introduced as heat sources in the thermal static solver. The copper temperature is used to compute the coils resistors introduced in the magnetic solver and this loop continues until change in temperatures is less than 0.1 °C. Both solvers use the same mesh and are included in Opera3D software.

B. Analytical Model

In order to motivate the need of a time expensive 3D FEA model, an Analytical Model (AM) is also used to compare RBDO results with both models.

The physical phenomena within the transformer are electric, magnetic and thermal. The assumptions for AM are uniform distribution of the magnetic flux density in the iron core and no voltage drop due to the magnetizing current. The thermal assumptions are the same than the 3D FEA except that the temperatures are uniform within the coils and laminations.

The weakest points of AM are the assumption of uniform temperature in copper and iron, and the approximation of the leakage inductance values.

C. Optimization Problem

The optimization problem contains 7 design variables. There are three parameters (𝑎, 𝑏, 𝑐) for the shape of the lamination, one for the frame (𝑑), two for the section of conductors (𝑆 1 , 𝑆 2 ), and one for the number of primary turn (𝑛 1 ) (Fig. 3

, left).

There are 7 inequality constraints in this problem. The copper and iron temperatures 𝑇 𝑐𝑜 , 𝑇 𝑖𝑟 respectively should be less than 120℃ and 100℃. The efficiency 𝜂 should be greater than 80% . The magnetizing current 𝐼 𝜇 𝐼 1 ⁄ and drop voltage Δ𝑉 2 𝑉 2 ⁄ should be less than 10% . All these constraints are computed with FEM or AM model. Finally, the filling factors of both coils 𝑓 1 , 𝑓 2 should be lower than 0.5.

The goal is to minimize the mass 𝑚 𝑡𝑜𝑡 of iron and copper materials. Thus, the optimization problem is expressed as: min 𝑚 𝑡𝑜𝑡 (𝑎, 𝑏, 𝑐, 𝑑, 𝑆 1 , 𝑆 For RBDO, all constraints are considered with a target probability of failure equal to 0.13%, which means a reliability index of 3. The standard deviation of each design variable is equal to 1% of its lowest bound.

D. Results

Table 2 shows optimal values, objective, probabilities of failure calculated by MCS with 10 6 samples computed with the meta-model, and the number of evaluations. For AM, SORA without meta-model is also tested and the number of evaluations is greater than 10,000, so it can be seen that SORA3 with kriging meta-model (SORA3 + AM, center column in Table 2) can find a solution almost satisfying all constraints with less evaluations. However when the same solution is reevaluated with FEM (FEM reeval., right column in Table 2), the highest probability of failure is 90%, so RBDO cannot be performed with AM only. The mass computed with FEM is slightly different from the one with AM because the voltage drop is considered to calculate the number of turns for the secondary coil.

SORA3 with FEM (SORA3 + FEM, left column in Table 2) leads to a probability of failure close to its target value. The objective value is higher with FEM because AM underestimates constraints. The initial sampling includes 7,000 points evaluated in parallel on 24 cores in about 49 hours then the 265 infill sampling points are evaluated sequentially in about 44 hours. The first advantage of SORA3 with FEM is that a significant computing time can be saved as it reduces the number of evaluations. The second advantage is that kriging model gives accurate derivatives that enable the use of fast gradient-based algorithm. Contrarily, as FEM provides noisy derivatives it requires a noise-free costly algorithm when directly connected with it. VI. CONCLUSION According to the mathematical example, the third strategy of kriging-based SORA is the most efficient without losing too much accuracy among the 6 approaches proposed.

RBDO of a single-phase safety isolating transformer is also performed here with FEM and the kriging-based SORA shows its applicability in dealing with this highly constrained problem by reducing the number of evaluations. Then, compared with analytical model of the same device, this approach with FEM could get a more accurate solution.

Fig. 1 .

 1 Fig. 1. The process of SORA3 strategy.

Fig. 3 .

 3 Fig. 3. Design variables of transformer and mesh.

Fig. 2 .

 2 Fig. 2. Iterations of SORA3 for the mathematical example (green points are initial sampling, pink points are enrichment samples during deterministic optimizations, yellow ones are added by reliability analysis and blue ones are MPTPs at the current iteration, dashed lines and contours present the real constraints and objectives respectively while solid ones present the metamodels).

TABLE 1 RESULTS

 1 OF MATHEMATICAL EXAMPLE USING DIFFERENT STRATEGIES

	Strategy	Number of evaluations	Optimal solution	Optimal value	Maximal 𝑃 𝑓 (%)
	SLA (exact model)	165	[2.2512; 1.9677] -1.9953 2.32
	PMA/SORA (exact model)	3183/531	[2.2513; 1.9691] -1.9945 2.27
	SLA	26	[2.2466; 1.9617] -1.9996 2.59
	PMA1	29	[2.2494; 1.9649] -1.9972 2.44
	PMA2	1804	[2.2513; 1.9691] -1.9945 2.27
	SORA1/2/3	142/97/45	[2.2513; 1.9691] -1.9945 2.27

TABLE 2

 2

	RESULTS OF TRANSFORMER OPTIMIZATION WITH META-MODEL
	Values	SORA3 + FEM SORA3 + AM FEM reeval.
	𝑎	12.902	13.153
	𝑏	46.042	51.039
	𝑐	18.183	16.532
	𝑑	42.318	43.098
	𝑛 1	659.06	641.75
	𝑆 1	0.3254	0.3216
	𝑆 2	2.7552	2.8956
	𝑚 𝑡𝑜𝑡	2.4028	2.3552	2.3520
	𝑃(𝑇 𝑐𝑜 > 120℃)	0%	0%	0%
	𝑃(𝑇 𝑖𝑟 > 100℃) ⁄ > 0.1) 𝑃(∆𝑉 2 𝑉 2 𝑃(𝐼 μ 𝐼 1 ⁄ > 0.1)	0.1506% 0% 0.1348%	0.1567% 0% 0.1420%	90.05% 0% 0.3281%
	𝑃(𝑓 1 > 0.5)	0.1327%	0.1236%	71.20%
	𝑃(𝑓 2 > 0.5)	0.1282%	0.1307%	0.0014%
	𝑃(𝜂 < 0.8)	0%	0%	0%
	Evaluations	7265	7242	1


Professor Harry Sokol
email: harry.sokol@aphp.fr

MD Nizar Mahlaoui

MD, PhD Claire Aguilar

MSc Perrine Bach

MD, PhD Olivier Join-Lambert

Aurélie Garraffo

MD, PhD Philippe Seksik

François Danion

MSc Sarah Jegou

MSc Marjolene Straube

MD François Danion

MSc Christelle Lenoir

MD, PhD Bénédicte Neven

Despina Moshous

PhD Stéphane Blanche

MD, PhD Bénédicte Pigneur

MD, PhD Olivier Goulet

MD, PhD Frank Ruemmele

MD, PhD Felipe Suarez

MD, PhD Laurent Beaugerie

MD, PhD Stéphanie Pannier

MD, PhD Françoise Mazingue

MD, PhD Olivier Lortholary

MD, PhD Lionel Galicier

MD, PhD Capucine Picard

MD, PhD Geneviève De Saint Basile

PhD Sylvain Latour

Professor Alain Fischer
email: alain.fischer@aphp.fr

Galicier

Francoise Mazingue

François Danion

Francois Danion

Intestinal dysbiosis in Inflammatory Bowel Disease associated with primary immunodeficiency

Keywords: gut microbiota, inflammatory bowel disease, primary immunodeficiency

To the Editor:

The gut microbiota is composed of a huge amount and diversity of microorganisms playing major roles in physiological and pathological settings such as inflammatory bowel disease (IBD). The triggering role of the gut microbiota on intestinal inflammation is known and an intestinal dysbiosis (i.e. an imbalance in the microbiota composition and function) has been pointed out in IBD 15 years ago. Moreover, a positive efficacy signal has been observed in studies evaluating fecal microbiota transplantation in IBD demonstrating that the gut microbiota alterations is an actor in the inflammatory process and not simple consequence [START_REF] Pigneur | Fecal microbiota transplantation in inflammatory bowel disease: the quest for the holy grail[END_REF] .

Environmental factors have an important effect on gut microbiota composition, but host genetic has an impact too. This effect has been demonstrated in knockout mouse for innate immunity genes. Microbiome genome-wide association studies have identified several associations between genetic polymorphisms and the gut microbiota composition offering some clues for the effects of genes on microbiota composition in humans [START_REF] Hall | Human genetic variation and the gut microbiome in disease[END_REF] . However, the effect of major specific disease-causing genetic defects on gut microbiota has never been studied in humans yet.

Here we characterized the fecal microbiota composition of patients with three types of rare primary immunodeficiency (PID) causing IBD conditions, chronic granulomatous disease (CGD, 11 samples), X-linked inhibitor of apoptosis (XIAP, 7 samples) deficiency and partial Tetratricopeptide Repeat Domain 7A (TTC7A, 7 samples) deficiency, in comparison to patients with non-genetic determined IBD (18 samples) and healthy subjects (HS, 23 samples, Suppl Table 1).

The microbiota composition was assessed by 16S sequencing. As seen previously, gut microbiota composition was influenced by age, antifungal and antibiotic treatment 3 (Suppl Figure 1). Beta diversity analysis showed a remarkable clustering of samples according to disease phenotypes (Figure 1A). The difference with HS was statistically significant for CGD and TTC7A groups while it did not reach significance for XIAP (Figure 1 A-D). Interestingly the difference with IBD patients was statistically significant for all three PID patients groups (Figure 1A). Age was very heterogeneous in patients from TTC7A and XIAP groups (0.62 to 18 and 1.3 to 34 years respectively). The gut microbiota composition changes constantly between birth and the age of 3 years when it reaches an adult-like configuration [START_REF] Tamburini | The microbiome in early life: implications for health outcomes[END_REF] . We thus performed the same analysis after segregating the patients according to their age (younger or older than 3 years). The difference with HS was maintained for TTC7A patients in both age groups and reached significance in XIAP patients older than 3 (only 2 XIAP patients younger than 3) (Suppl Figure 2).

Patients with CGD, XIAP and TTC7A deficiency can experience intestinal inflammation mimicking IBD [START_REF] Uhlig | Monogenic diseases associated with intestinal inflammation: implications for the understanding of inflammatory bowel disease[END_REF] . In our study population, samples were taken from patients with (PIDIBD) or without (PIDno-IBD) IBD involvement. In addition, PID patients could be in active phase (flare) or remission of IBD (PID, IBD flare/inactive). Microbiota from patients with active intestinal inflammation (PIDIBD (flare)) was significantly different from the one of PIDno-IBD patients (Figure 1E). We observed a significant diversity reduction in patients with PIDIBD compared to HS (Figure 1F, Suppl Figure 3A). Moreover, this reduction was even stronger in patients with active intestinal inflammation (Figure 1G, Suppl Figure 3B). These effects were seen in the three-studied PID although the numbers were too low to reach statistical significance when analyzing each PID independently (Suppl Figure 4).

In accordance with the published literature, the bacterial microbiota of all groups was dominated by bacteria from Firmicutes, Bacteroidetes, Actinobacteria and Proteobacteria phyla. Differences were noticed between PID patients and HS with notably a dramatic increase in Proteobacteria from the Enterobacteriaceae family in the TTC7A group and an increase proportion in bacteria from the Bacteroidetes phylum and from the Clostridiaceae family in the CGD and XIAP groups respectively (Figure 1 H, Suppl Figure 3C). These alterations were found both in patients younger and older than 3 years (Suppl Figure 4 B,D) and were clearly stronger in patients with IBD and PIDIBD with active inflammation (Figure 1I; Suppl Figure 4 C, E). We then looked for the microbial features associated with the three studied PID and used a multivariate association test (MaAsLin) to control for potential confounding factors such as age, gender, smoking, intestinal inflammation and treatment. We observed several alterations specific to each PID compared to HS. Patients with TTC7A deficiency exhibited an increased abundance of Proteobacteria involving several components of the Gammaproteobacteria and Epsilonproteobacteria classes (Figure 2A). On the other hand, bacteria from the Ruminococcaceae family and notably the Oscillospira genera were decreased. Patients with XIAP deficiency exhibited an increased abundance of several bacterial taxa including members of the Proteobacteria, Firmicutes, Actinobacteria and Fusobacteria phyla (Figure 2B). Four of these taxa (Scardovia, Fusobacterium, Rothia dentocariosa, Veillonella) are not usually found in the gut but are known members of the oral microbiota that are implicated in dental caries, but also IBD, colorectal cancer, and liver diseases [START_REF] Kummen | The gut microbial profile in patients with primary sclerosing cholangitis is distinct from patients with ulcerative colitis without biliary disease and healthy controls[END_REF] . Interestingly, an increased abundance of oral bacteria in the gut has been observed in several intestinal and extra-intestinal diseases [START_REF] Atarashi | Ectopic colonization of oral bacteria in the intestine drives TH1 cell induction and inflammation[END_REF] .

Moreover, a recent study showed that some members of the oral microbiota can have proinflammatory effects when colonizing the gut [START_REF] Atarashi | Ectopic colonization of oral bacteria in the intestine drives TH1 cell induction and inflammation[END_REF] suggesting a possible direct effect of these oral bacteria in the intestinal inflammatory phenotype observed in patients with XIAP deficiency.

Another important specificity in XIAP patients was the presence of Lactococcus garvieae which is a highly virulent pathogen affecting saltwater fish [START_REF] Vendrell | Lactococcus garvieae in fish: a review[END_REF] and is rarely involved in human infection. L. garvieae was found at a high level (0.23-0.5% of all reads) in XIAP patients with active intestinal inflammation but not in all other patients in this study. Patients with CGD exhibited an increased abundance of Ruminococcus gnavus that has also been associated with ileal Crohn's disease [START_REF] Sokol | Fungal microbiota dysbiosis in IBD[END_REF] (Figure 2C).

In summary, we showed that the gut microbiota of patients with CGD, TTC7A and XIAP deficiency have distinct alterations suggesting a primary defect in host immune system as a basis of dysbiosis. Although it remains to be experimentally documented, the microbial alterations induced by the host genetic defect might play a role in some aspects of the PID phenotype and particularly intestinal involvement. The number of subjects studied was low because of the extreme rarity of the studied diseases. However, independent and ideally larger and longitudinal studies are required to confirm our findings. Controlling for potential confounding factors such as diet, treatment, age, delivery mode, socio-economic features, gastrointestinal symptoms, is particularly difficult in studies aiming at discriminating genetic factors from environmental and inflammation influences. Similarly, identification of ideal controls is not trivial. Although non-affected healthy siblings sharing similar lifestyles would be attractive, there are obvious limitations in term of feasibility.

Finally, if these observations are confirmed, the alteration in gut microbiota composition might have clinical interest as diagnosis biomarkers.

Harry Sokol 1,2,3,4 , Nizar Mahlaoui 4,5,17 , Claire Aguilar 6,8,17 , Perrine Bach 4 , Olivier Join-Lambert 7 , Aurélie Garraffo ,6 , Philippe Seksik 1,3 , François Danion 6 , Sarah Jegou 1 , Marjolene Straube 2 , Christelle Lenoir 8 , Bénédicte Neven 5,8,9 , Despina Moshous 5,8,9 , Stéphane Blanche 5,9 , Bénédicte Pigneur 8,10 , Olivier Goulet 8,10 , Frank Ruemmele 8,10 , Felipe Suarez [START_REF] Tamburini | The microbiome in early life: implications for health outcomes[END_REF][START_REF] Vendrell | Lactococcus garvieae in fish: a review[END_REF]12 , Laurent Beaugerie 3 , Stéphanie Pannier 11,13 , Françoise Mazingue 14 , Olivier Lortholary [START_REF] Tamburini | The microbiome in early life: implications for health outcomes[END_REF][START_REF] Kummen | The gut microbial profile in patients with primary sclerosing cholangitis is distinct from patients with ulcerative colitis without biliary disease and healthy controls[END_REF][START_REF] Sokol | Fungal microbiota dysbiosis in IBD[END_REF] , Lionel

Figure 1 :

 1 Figure 1: Abnormal bacterial microbiota composition and diversity in PID. Principal

Figure 2 :F

 2 Figure 2: Bacterial taxa associated with TTC7A, XIAP and CGD. Bacterial taxa associated

  n o -I B D P I D I B D ( n o n -a c t i v e ) P I D I B D ( f l a r e )

Age: median (min-max) 9 (0.83-32) 28.5 (10-41) 10 (6-31) 1.8 (0.62-18) 8 (1.3-34) Male Gender % ( n) 61% ( 14) 89% ( 16) 100% (11) 43% (3) 100% (7) Gut inflammation history 0 100% (18) 82% ( 9) 100% ( 7) 43%

Active gut inflammation 0 0 64% (7) 71% ( 5

SUPPLEMENTARY INFORMATION

MATERIAL AND METHODS

Patients and sample collection

Patients with primary immunodeficiency were recruited at Necker-Enfants malades University IV, Suivitheque study). None of the IBD patients and the healthy subjects reported having taken antibiotics, probiotics or using colon-cleansing products for at least 1 months prior to enrollment. Patient characteristics are presented in Table 1. Whole stools were collected in sterile boxes and immediately homogenized, and 0.2 g aliquots were frozen at -80°C for further analysis. Patients' data were collected, including clinical, treatment, immunological and genetic diagnosis. Enrolled patients with PID had damaging causal mutation respectively in either the TTC7A gene "ELA syndrome" 26 , BIRC4 (XIAP deficiency) and CYBBNox2 (X-linked CGD).

None of the subjects received probiotics.

Fecal DNA extraction

Genomic DNA was extracted from 200 mg of feces as previously described. 27 . Following microbial lysis with both mechanical and chemical steps, nucleic acids were precipitated in isopropanol for 10 minutes at room temperature, incubated for 15 minutes on ice and centrifuged for 30 minutes at 15,000 g and 4°C. Pellets were suspended in 112 µL of phosphate buffer and 12 µL of potassium acetate. After RNase treatment and DNA precipitation, nucleic acids were recovered via centrifugation at 15,000 g and 4°C for 30 minutes. The DNA pellet was suspended in 100 µL of TE buffer.

16S rRNA gene sequencing

Fecal DNA was extracted from the weighted feces before and during the infection as previously described 27 . Microbial diversity was determined for each sample by targeting a portion of the ribosomal genes. A 16S rRNA gene fragment comprising V3 and V4 hypervariable regions (16S; 5′-TACGGRAGGCAGCAG-3′ and 5′-CTACCNGGGTATCTAAT-3′) was amplified using an optimized and standardized 16S-amplicon-library preparation protocol (Metabiote, GenoScreen). Briefly, 16S rRNA gene PCR was performed using 5 ng genomic DNA according to the manufacturer's protocol (Metabiote) using 192 bar-coded primers (Metabiote MiSeq Primers, GenoScreen) at final concentrations of 0.2 μM and an annealing temperature of 50°C for 30 cycles. The PCR products were purified using an Agencourt AMPure XP-PCR Purification system (Beckman Coulter), quantified according to the manufacturer's protocol, and multiplexed at equal concentrations. Sequencing was performed using a 250-bp paired-end sequencing protocol on an Illumina MiSeq platform (Illumina) at GenoScreen. Raw paired-end reads were subjected to the following process: (1) quality filtering using the PRINSEQ-lite PERL script 5 by truncating the bases from the 3′ end that did not exhibit a quality < 30 based on the Phred algorithm; (2) paired-end read assembly using FLASH (fast length adjustment of short reads to improve genome assemblies) with a minimum overlap of 30 bases and a 97% overlap identity; and (3) searching and removing both forward and reverse primer sequences using CutAdapt, with no mismatches allowed in the primers sequences. Assembled sequences for which perfect forward and reverse primers were not found were eliminated.

16S rRNA gene sequence analysis

The sequences were demultiplexed and quality filtered using the Quantitative Insights Into Microbial Ecology (QIIME, version 1.9.1) software package 28 , and the forward and reverse Illumina reads were joined using the fastq-join method (http://code.google.com/p/ea-utils). The sequences were assigned to OTUs using the UCLUST algorithm 29 with a 97% threshold of pairwise identity and classified taxonomically using the Greengenes reference database 30 .

Principal component analyses of the Bray Curtis distance were built and used to assess the variation between experimental groups (beta diversity). Significance was assessed using ANOSIM (9999 permutations). The Shannon and Chao1 diversity index were calculated using rarefied data (depth = 17,000 sequences/sample) and used to characterize species diversity in a community.

Statistical analysis

GraphPad Prism version 6.0 (San Diego, CA) was used for all analyses and graph preparation.

For all graph data, the results are expressed as the mean ± SEM, and statistical analyses were performed using the 2-tailed nonparametric Mann-Whitney U-test or Kruskal-Wallis test with Dunn's Multiple Comparison Test. Statistical significance of sample grouping for beta diversity analysis was performed using Anosim method (9999 permutations). Differences with a p value less than 0.05 were considered significant.

Multivariate Analysis by Linear Models (MaAsLin), a multivariate statistical framework, was used to find associations between clinical metadata and microbial community abundance 31 . We took into account as much factor as possible in the analysis to control for the effects of potential confounding factors including: age, gender, smoking, intestinal inflammation and treatment and enteral nutrition.

SUPPLEMENTARY FIGURE LEGEND


François Charles
email: francois.charles@obs-banyuls.fr

Pierre-Guy Sauriau

Fabien Aubert

Benoît Lebreton

François Lantoine

Pascal Riera

Sources partitioning in the diet of the shipworm Bankia carinata (J.E. Gray, 1827): an experimental study based on stable isotopes

Keywords: Teredinidae, Trophic ecology, Stable isotope, in situ experiment, Mediterranean Sea Graphical

Adaptations that allow teredinids to maintain and thrive on wood, a nutritionally unbalanced food, make these marine bivalves remarkable. Capable of filter-feeding, shipworms house endosymbiotic bacteria synthesizing cellulolytic enzymes for digestion of wood carbohydrates and providing nitrogen to their host through nitrogen fixation. To what extent each of these nutrition modes contributes to the shipworm's metabolism remains an open question. In this experimental study, we estimated source partitioning through the determination of δ 13 C and δ 15 N values in original biological samples. For this purpose, pieces of common alder (Alnus glutinosa) were immersed at a coastal station of the northwestern Mediterranean Sea. The shipworm Bankia carinata infected wood logs and stable isotope mixing models suggested it got most of the carbon and nitrogen it needs from separate sources. From 71 to 77% of the carbon was derived from the digestion of wood carbohydrates, whereas between 42 and 82% of the nitrogen originated from N2 fixation. These first semi-quantitative estimations suggest that the contribution of N2 fixers to nitrogen requirements of this shipworm species is far from incidental.

Introduction

Wood-boring teredinid bivalves occur on wood scattered at sea. In spite of being pests for all kinds of maritime wooden structures, these mollusks, commonly known as shipworms, are one of the finest examples of adaptation within the tree of life [START_REF] Cragg | Lignocellulose degradation mechanisms across the Tree of Life[END_REF]. The strategy that enables these bivalves to digest wood and supplement a nutritionally imbalanced diet has been debated for over a century. While it was shown some time ago that shipworms were able to digest carbohydrates from wood [START_REF] Dore | The digestion of wood by Teredo navalis[END_REF][START_REF] Miller | Digestion of wood by shipworms[END_REF], the coverage of nitrogen requirements has remained under discussion [START_REF] Carpenter | Nitrogen fixation in marine shipworms[END_REF][START_REF] Greenfield | Observation on the nitrogen and glycogen content of Teredo (Lyrodus) pedicellate de Quatrefages at Miami, Florida[END_REF][START_REF] Lasker | The origin and distribution of nitrogen in Teredo bartschi Clapp[END_REF][START_REF] Nishimoto | Organic carbon and nitrogen source of sunken wood communities on continental shelves around Japan inferred from stable isotope ratios[END_REF][START_REF] Pechenik | The role of phytoplankton in the diets of adult and larval shipworms, Lyrodus pedicellatus (Bivalvia: Teredinidae)[END_REF][START_REF] Yamanaka | Re-evaluation of nutrient sources for deep-sea wood-boring bivalves using the isotopic composition of bulk C, N, S, and amino acid nitrogen[END_REF] to explain how some species [START_REF] Gallager | Physiological aspects of wood consumption, growth, and reproduction in the shipworm Lyrodus pedicellatus Quatrefages (Bivalvia: Teredinidae)[END_REF] possess the ability to maintain and thrive in the absence of any other food source than wood.

Rather than taking advantage of an association with its gut microflora which otherwise is scarce [START_REF] Betcher | Microbial distribution and abundance in the digestive system of five shipworm species (Bivalvia: Teredinidae)[END_REF], shipworms produce endogenous lignocellulolytic enzymes [START_REF] Sabbadin | Uncovering the molecular mechanisms of lignocellulose digestion in shipworms[END_REF], but also rely on dense populations of intracellular bacterial symbionts hosted in the gland of Deshayes, a cluster of bacteriocytes located within the gills [START_REF] Popham | Bacterial associations in the Teredo Bankia australis (Lamellibranchia: Mollusca)[END_REF]. The many types of bacteria housed in these eukaryotic cells can both digest wood carbohydrates and fix nitrogen gas [START_REF] Distel | Teredinibacter turnerae gen. nov., sp. nov., a dinitrogen-fixing, cellulolytic, endosymbiotic gamma-proteobacterium isolated from the gills of wood-boring molluscs (Bivalvia: Teredinidae)[END_REF][START_REF] Waterbury | A cellulolytic nitrogen-fixing bacterium cultured from the gland of Deshayes in shipworms (Bivalvia: Teredinidae)[END_REF]. It has been demonstrated that wood-degrading enzymes produced by the endosymbionts were selectively translocated from the gills to the gut of the host (O' Connor et al., 2014), and that dissolved nitrogen fixed by the intracellular bacterial symbionts within the bacteriocytes could further be used for host metabolism [START_REF] Lechene | Quantitative imaging of nitrogen fixation by individual bacteria within animal cells[END_REF]. All this makes teredinids a unique model in animal endosymbiosis [START_REF] Lechene | Quantitative imaging of nitrogen fixation by individual bacteria within animal cells[END_REF] and may explain how shipworms are able to use wood as a primary food source.

Marine wood-boring bivalves are indeed remarkable consumers in that the number of potential feeding sources differs according to whether one considers their carbon or nitrogen requirement. With regard to carbon sources, teredinids can rely on wood but also, as filter-feeding bivalves, on suspended particulate organic matter from sea water. Regarding nitrogen sources, wood content is seemingly inadequate to meet shipworms' requirements [START_REF] Lasker | The origin and distribution of nitrogen in Teredo bartschi Clapp[END_REF][START_REF] Carpenter | Nitrogen fixation in marine shipworms[END_REF], but nitrogen incorporation from this source cannot be excluded [START_REF] Greenfield | Observation on the nitrogen and glycogen content of Teredo (Lyrodus) pedicellate de Quatrefages at Miami, Florida[END_REF][START_REF] Yamanaka | Re-evaluation of nutrient sources for deep-sea wood-boring bivalves using the isotopic composition of bulk C, N, S, and amino acid nitrogen[END_REF]. Suspended organic particles are also sources of nitrogen but may remain incidental [START_REF] Gallager | Physiological aspects of wood consumption, growth, and reproduction in the shipworm Lyrodus pedicellatus Quatrefages (Bivalvia: Teredinidae)[END_REF] because the ability to filter feed in most shipworm species is reduced by the shortening of gill lamellae [START_REF] Turner | A survey and illustrated catalogue of the Teredinidae (Mollusca: Bivalvia)[END_REF][START_REF] Distel | Molecular phylogeny of Pholadoidea Lamarck, 1809 supports a single origin for xylotrophy (wood feeding) and xylotrophic bacterial endosymbiosis in Bivalvia[END_REF]. Nitrogen fixation by endosymbionts [START_REF] Lechene | Quantitative imaging of nitrogen fixation by individual bacteria within animal cells[END_REF] is one further source. Even though the routes of C and N are known and common among shipworm species, the relative contribution of these routes to the shipworms' metabolism remains controversial particularly with respect to nitrogen.

Stable isotope compositions of marine organic matter, terrestrial plants and N2 dissolved in seawater are different [START_REF] Peterson | Stable isotopes in ecosystem studies[END_REF], and can therefore help in deciphering the composition of the shipworm's diet. To our knowledge, the few studies published on this topic were not conclusive about source partitioning [START_REF] Nishimoto | Organic carbon and nitrogen source of sunken wood communities on continental shelves around Japan inferred from stable isotope ratios[END_REF]; Paalvast and van der Velde, 2013; [START_REF] Yamanaka | Re-evaluation of nutrient sources for deep-sea wood-boring bivalves using the isotopic composition of bulk C, N, S, and amino acid nitrogen[END_REF]. In this context, the aim of the present study was thus to investigate carbon and nitrogen source partitioning by shipworms via an in situ experimental approach and through the use of stable isotopes (δ 13 C, δ 15 N).

Material and methods

Experimental design

The experiment was performed from June 6 to October 24, 2016 in the Bay of Banyuls-sur-Mer, in the north-western Mediterranean Sea, France. Twelve 15 cm long sections were sawn from, a branch of a common alder tree, Alnus glutinosa (L.) Gaertner, a widespread tree species along the banks of La Massane, a local river which drains a small forested basin on the foothills of the French Pyrenees mountains. The bark-intact small logs were air-dried indoors for a full year before being arranged into packs and deployed at a coastal station (26 m depth; 42°29'302 N; 03°08'700 E). Each pack consisted of four 2.5 to 5.5 cm cross-section wooden pieces tied together with a Colson™ collar. Two packs were installed in June and an additional one on August 20. Wood packs were attached at one end to an anchor and at the other end to a polyethylene buoy in such a way that the wood blocks remained suspended between 1 and 1.5 m above the seabed (Fig. S1).

Sampling and preparation for stable isotope analyses

Potential food sources, namely alder wood within the experimental wood packs and suspended particulate organic matter (SPOM) from the surrounding sea water, were sampled on several occasions during the experiment. Alder wood was sampled on June 6, at the time of deployment of the first wood packs, and on October 24 after time intervals of nearly 2 and 4.5 months from submersion. SPOM was sampled when installing the experiment, at two intermediate dates, July 7 th and September 8 th , and at the end of the experiment. Consumers, namely animals present either on (i.e. Mytilus galloprovincialis, Hiatella arctica, Musculus subpictus, Serpula vermicularis, Spirobranchus triqueter, Ascidia conchilega) or inside (Bankia carinata) the wood logs, were sampled once, on October 24.

Three water samples were collected for SPOM near the submerged wood with a Niskin bottle (10 L) and pre-filtered through a 200 μm mesh to remove large zooplankton and detritus. SPOM was then recovered by filtration on pre-combusted Whatman GF/F filters, 25 mm in diameter, 0.7 µm pore size. The filters were then placed in individual petri dishes and dried in an oven at 50°C overnight. Dry filters were kept away from light in a desiccator until analysis.

Invertebrate species present at the surface of the wood packs were collected manually. Wood logs were then cut with a chisel for the collection of wood samples and shipworms. All the animals present were kept alive for 24h in filtered sea water to allow evacuation of the gut contents and, if necessary, were then dissected to separate flesh from the shells. Wood chips and animal tissues were rinsed with deionised water and stored in aluminium foil at -20°C until analysis. Freeze-dried tissues and wood samples were ground into a fine powder using a ball mill (MM 400, Retsch) and a sander (R3000, Dremel), respectively. These samples were then weighed into tin capsules (ca. 0.3 to 0.4 mg for all sample tissue and wood carbon analysis, and ca. 4 to 5 mg for wood nitrogen analysis).

SPOM samples may contain carbonates that are more enriched in 13 C than plant and animal tissues [START_REF] Deniro | Influence of diet on the distribution of carbon isotopes in animals[END_REF]. Accordingly, δ 13 C and δ 15 N were analysed separately. For measurements of δ 13 C values, one half of each filter was acidified for 4 h with 37 % HCl acid vapours under moderate vacuum in a glass desiccator to remove carbonates (Malet et al., 2007). For measurements of δ 15 N values, the remaining half of each filter was analysed without any pre-treatment, as acidification may alter nitrogen isotope composition [START_REF] Bunn | Effects of acid washing on stable isotope ratios of C and N in penaeid shrimp and seagrass: implications for food-web studies using multiple stable isotopes[END_REF][START_REF] Kennedy | The effect of acidification on the determination of organic carbon, total nitrogen and their stable isotopic composition in algae and marine sediment[END_REF]. The two sets of sub-samples (i.e. acidified half-filters and non-acidified half-filters) were packed separately into individual tin capsules.

Stable isotope ratio measurements

Carbon and nitrogen isotope composition was determined using an elemental analyser (Flash EA 1112, Thermo Scientific, Milan, Italy) coupled to a continuous-flow isotope-ratio mass spectrometer (Delta V Advantage with a Conflo IV interface, Thermo Scientific, Bremen Germany). Analyses were conducted at the LIENSs stable isotope facility at the University of La Rochelle, France. Data are expressed in the δ notation (in ‰) as deviations from international standards (Vienna Pee Dee Belemnite for δ 13 C and N2 in air for δ 15 N) following the formula: δX = [(Rsample /Rreference ) -1] × 10 3 ,with R = 13 C/ 12 C for carbon and 15 N/ 14 N for nitrogen. Calibration was done using reference materials (USGS-24, IAEA-CH6, IAEA-600 for carbon; IAEA-N2, IAEA-NO-3, IAEA-600 for nitrogen). Analytical precision was 0.15 ‰ for both C and N isotope analyses based on the analyses of acetanilide (Thermo Scientific) and peptone (Sigma-Aldrich) used as internal laboratory standards.

Data analysis and statistical treatment

All statistics were performed using R, released by the R Foundation for Statistical Computing (R Core Team, 2013). First, from the R commander package (Fox and Vialat, 2017), non-parametric Mann-Whitney U tests were run for testing the hypothesis that the distribution of data was the same between groups. In this context, data accounted for either δ values or C/N ratio; groups for exposure intervals to seawater (2 vs 4.5 months), food sources (wood vs SPOM), or feeding modes (suspensionfeeding vs xylotrophy).

The Bayesian mixing model SIAR (Stable Isotope Analysis in R; Parnell et al., 2010) was used to infer the feasible contributions of carbon and nitrogen sources to the diet of the shipworms. Models were run for 200,000 iterations and the first 50,000 iterations were discarded. For SIAR calculations, A. glutinosa wood and SPOM were considered as potential sources of carbon while A. glutinosa, SPOM and dissolved N2 were considered as potential sources of nitrogen. Isotope compositions of SPOM and wood used as end-members into mixing models were averages of values measured all along the study to account for potential temporal variations of their isotope compositions. Trophic fractionation factors used in mixing models were 0.30 ± 0.21‰ for δ 13 C values and 2.5 ± 0.25‰ for δ 15 Concerning the specific utilization of dissolved N2 by the bivalves, we thus used the δ 15 N of 0.6‰ given by [START_REF] Sigman | Nitrogen isotopes in the ocean[END_REF] for the marine dissolved N2 and a fractionation of -1‰ as suggested by these authors for the fixed N input to the ocean from N2 fixation. No fractionation was expected between the endosymbionts and their host [START_REF] Conway | The role of endosymbiotic bacteria in the nutrition of Solemya velum: Evidence from a stable isotope analysis of endosymbionts and host[END_REF]).

Results

Stable isotope compositions

δ 15 N and δ 13 C values of alder wood remained in narrow ranges around mean values of -2.2‰ and -29.2‰ for N and C isotopes, respectively (Table 1, Figure 1). Ranges of SPOM isotope compositions were slightly wider around mean time-integrated values of 3.4‰ and -23.8‰ for δ 15 N and δ 13 C values, respectively. Alder wood was significantly more 15 N depleted and 13 C depleted than SPOM and contained two orders less nitrogen relative to its carbon content compared to SPOM (one-sided Mann-Whitney U-test, p = 1 in all cases).

The number of animals present either on or inside the wood logs increased with exposure time intervals. Only one species of shipworm (Table 2), namely, Bankia carinata (J.E. Gray 1827), was observed. Most other taxonomic groups consisted of suspension-feeding species. As no significant difference for both δ 15 N and δ 13 C values was observed between the two exposure intervals to seawater (two-sided Mann-Whitney U test, p = 0.09, p = 0.10), B. carinata was considered a homogenous group.

The average δ 15 N value of the shipworm bodies was close to 0‰ and their content was highly depleted in 13 C with a mean δ 13 C value of -27.5‰. Shipworms were much more 15 N and 13 C depleted than known strict seston feeders such as M. galloprovincialis, H. arctica, M. subpictus, S. vermicularis, S. triqueter or A. conchilega (one-sided Mann-Whitney U-tests, p = 1 for all cases). Higher C/N ratios were observed in shipworms (one-sided Mann-Whitney U-tests, p = 1) due to higher carbon (39.9 ± 1.9% versus 34.9 ± 3.3%) and lower nitrogen contents (5.2 ± 1.6% versus 8.4 ± 1.6%) than in other consumers. The addition of trophic fractionation factors to alder wood isotope compositions demonstrates a slight shift between expected and measured of δ 13 C values in shipworms (Figure 1).

Stable isotope mixing models

Based on 95% Bayesian credible intervals, wood and SPOM proportions to shipworms carbon diet ranged from 71 to 77% (mode = 74%), and 23 to 29% (mode = 26%), respectively (Figure 2A). As a result, alder wood contributed to most of the carbon assimilated by the shipworms during the experiment as compared to local SPOM. Among the three nitrogen sources (Figure 2B), the Bayesian mixing model suggested a major contribution from N2 fixation (ranging from 42 to 82%, mode = 54%) and from A. glutinosa (from 9% to 51%, mode = 37%), much higher than contribution of SPOM (7 to 13%, mode = 9%).

Discussion

Bankia carinata resorted to all the pathways of carbon and nitrogen acquisition described in shipworms. The approach permitted us to determine estimates of nitrogen source contribution, even though the number of sources was too large to get a unique solution from standard linear mixing models.

As expected, the conversion of isotope compositions into contributions to the shipworm diet confirmed that wood was the main source of carbon, but it primarily revealed the significant contribution of nitrogen fixation pathways to the nitrogen content.

Source partitioning in the diet of B. carinata

Bankia carinata possesses a large caecum storing wood particles and all features associated with xylotrophy in shipworms, including the presence of dense populations of endosymbiotic bacteria within the large interlamellar space of the single demibranch gills [START_REF] Distel | Molecular phylogeny of Pholadoidea Lamarck, 1809 supports a single origin for xylotrophy (wood feeding) and xylotrophic bacterial endosymbiosis in Bivalvia[END_REF]. In European waters, this shipworm has only been reported in the Mediterranean Sea [START_REF] Turner | A survey and illustrated catalogue of the Teredinidae (Mollusca: Bivalvia)[END_REF] probably because of a preference for hypersaline conditions [START_REF] Borges | Diversity, environmental requirements, and biogeography of bivalve wood-borers (Teredinidae) in European coastal waters[END_REF]. Bankia carinata is a broadcast spawning species that releases fertilized eggs into the water for long periods of development, increasing larval dispersion. Growth in the genus Bankia is generally fast and produces adults of large size [START_REF] Haderlie | Settlement, growth rates and depth preference of the shipworm Bankia setacea (Tryon) in Monterey Bay[END_REF][START_REF] Macintosh | Contrasting life histories in shipworms: Growth, reproductive development and fecundity[END_REF]. This species is reported in the studied area together with Lyrodus pedicellatus, Teredo navalis and Nototeredo norvagica [START_REF] Charles | It's what's inside that counts: computer aided tomography for evaluating the rate and extent of wood consumption by shipworms[END_REF]. As population dynamics patterns depend on reproductive modes, differences must be expected in food source partitioning between shipworm species according to the level of exhaustion of the wood resource.

The measurement of isotope ratios on the whole body of the borers corresponded to the isotopic signature of the flesh of the shipworms including those of the bacteria from the digestive tract and the endosymbionts from the gills. The digestive tract, however, contains relatively few bacteria [START_REF] Betcher | Microbial distribution and abundance in the digestive system of five shipworm species (Bivalvia: Teredinidae)[END_REF] and endosymbiont bacteria represent a negligible amount of matter compared to the tissues of the shipworm. Thence, it is reasonable to assume that the measurements made on the shipworms were essentially the values of the isotopic ratios in the tissues of the host.

The differences of δ 13 C and δ 15 N values between B. carinata and SPOM highlighted that, despite being able to capture organic particles from water, B. carinata food resources drastically differ from those of strict filter feeders living in the same trophic conditions (Figure 1). Accordingly, from these results, most of the carbon in shipworms is clearly derived from the consumption of alder wood (Figure 2a). The carbon isotope mixing model suggested that SPOM could seemingly represent 25% of shipworm's diet. Structural tissues of wood mainly consist of carbohydrates and lignin but shipworms only digest carbohydrates [START_REF] Dore | The digestion of wood by Teredo navalis[END_REF][START_REF] Miller | Digestion of wood by shipworms[END_REF][START_REF] Sabbadin | Uncovering the molecular mechanisms of lignocellulose digestion in shipworms[END_REF].

Carbohydrates are more enriched in 13 C than lignin [START_REF] Benner | Depletion of 13 C in lignin and its implications for stable carbon isotope studies[END_REF], and this gives whole wood a δ 13 C value which probably led to underestimate the contribution of wood source. This was further suggested by the outputs of the mixing models based on δ 15 N values. High carbon content in the shipworms reflects the accumulation of large amounts of glycogen [START_REF] Greenfield | Observation on the nitrogen and glycogen content of Teredo (Lyrodus) pedicellate de Quatrefages at Miami, Florida[END_REF][START_REF] Potts | The structure and function of the liver of Teredo, the shipworms[END_REF] which serves as an energy store in bivalves (De Zwaan and Zandee, 1972) and has been shown to be exclusively derived from digestion of wood carbohydrates in the shipworm Lyrodus pedicellatus [START_REF] Lane | The distribution of glycogen content in the Shipworm Teredo (Lyrodus) pedicellata de Quatrefages at Miami, Florida[END_REF].

According to the SIAR mixing model, SPOM contributed about 10% in average of the nitrogen of the shipworms. Out of the two other potential sources, the mixing model suggested that nitrogen coming from endosymbiosis contributed to the nitrogen pool of the shipworm tissues at an equal or even greater level than alder wood. From the C/N ratio, it appears that B. carinata had reduced nitrogen requirements relative to other bivalves but a nitrogen content 20 times higher than that of alder wood.

Shipworms can rely on wood to acquire their nitrogen [START_REF] Lasker | The origin and distribution of nitrogen in Teredo bartschi Clapp[END_REF][START_REF] Yamanaka | Re-evaluation of nutrient sources for deep-sea wood-boring bivalves using the isotopic composition of bulk C, N, S, and amino acid nitrogen[END_REF] but also through the association with the N2-fixing endosymbionts hosted in their gills [START_REF] Lechene | Quantitative imaging of nitrogen fixation by individual bacteria within animal cells[END_REF]. Even though, by minimizing excretory losses, shipworms are extremely conservative with respect to nitrogen [START_REF] Gallager | Physiological aspects of wood consumption, growth, and reproduction in the shipworm Lyrodus pedicellatus Quatrefages (Bivalvia: Teredinidae)[END_REF], wood source alone can barely cover all their nitrogen needs. Together with recycling mechanisms, a polytrophic capacity is probably the best feeding strategy for balancing requirements in poor-nitrogen environments [START_REF] Muscatine | Reef corals: Mutualistic symbioses adapted to nutrient-poor environments[END_REF]. Intracellular bacteria populations from the gland of Deshayes produce cellulolytic enzymes which presumably in combination with enzymes produced endogenously [START_REF] Sabbadin | Uncovering the molecular mechanisms of lignocellulose digestion in shipworms[END_REF] allow shipworms to digest efficiently large amounts of wood carbohydrates. Endosymbiotic bacteria probably fulfill their N2 fixing function with an equal efficiency to the benefit of their host. Accordingly, the conversion of dissolved nitrogen gas into shipworm biomass is likely far from being incidental as suggested by the present study.

Carbon and nitrogen sources: scope of the study

The δ 15 N value of marine dissolved N2 was not measured but dissolved N2 in equilibrium with atmospheric N2 in the surface ocean has a δ 15 N value that does not vary greatly from 0.6‰ [START_REF] Sigman | Nitrogen isotopes in the ocean[END_REF]. Molecular nitrogen fixation commonly results in weak, slightly negative isotope fractionation (Wada and Hattori, 1979) and no fractionation between endosymbionts and their host was expected [START_REF] Conway | The role of endosymbiotic bacteria in the nutrition of Solemya velum: Evidence from a stable isotope analysis of endosymbionts and host[END_REF]. As a result, the δ 15 N value of dissolved gas nitrogen and fractionation factor implemented into the mixing model must be relatively close to the actual values. The slight 15 N enrichment of shipworm tissue relative to dinitrogen suggests that N2 fixation through endosymbiotic interaction contributed significantly to the nitrogen content of the shipworm. The contribution of other pathways of dinitrogen-fixing (i.g., ingestion of planktonic diazotrophs or wood particles from nitrogenfixing tree species) would imply trophic enrichment factors and that would have led to the higher 15 N enrichment of the shipworm.

Regarding wood source, Alnus glutinosa is a common tree species on the banks of European rivers. Alder wood, therefore, has good chance to end up at sea after events of river flooding. During this study, the characteristics of wood which was initially introduced dry into the water changed very little. δ 13 C values and C/N ratios were representative of most terrestrial C3 primary producers and in the range of those found in the literature for tree wood (Nishimoto et Zapata-Hernandez et al., 2016), including alder tree species [START_REF] Chambers | Nitrogen uptake and turnover in riparian woody vegetation[END_REF]. Nitrogen incorporation pathways differ among tree species. Alder, which grows on poor quality soils is, for instance, able to supplement its needs through an association with nitrogen-fixing bacteria hosted in its roots [START_REF] Hooker | The effectivity of Frankia for nodulation and nitrogen fixation in Alnus rubra and A. glutinosa[END_REF]. Alnus incana is known to derive almost all its nitrogen from the atmosphere [START_REF] Chambers | Nitrogen uptake and turnover in riparian woody vegetation[END_REF]. This assigned to the wood source a trophic position which probably led us to overestimate its contribution to the nitrogen content of the shipworm, at the expense of the other sources. It, therefore, should be kept in mind that the scope of this study concerns a nitrogenfixing tree species. δ 13 C and δ 15 N values of SPOM were consistent with previous measurements at the studied site [START_REF] Carlier | Benthic trophic network in the Bay of Banyuls-sur-Mer (northwest Mediterranean, France): An assessment based on stable carbon and nitrogen isotopes analysis[END_REF][START_REF] Nahon | Trophic relationships and UV-absorbing compounds in a Mediterranean medio-littoral rocky shore community[END_REF] with ranges of temporal changes of nearly 2‰ for both elements.

Marine SPOM is a mixture of living and dead particles [START_REF] Fry | Food web structure on Georges Bank from stable C, N, and S isotopic compositions[END_REF] of which the total amount and the individual contribution changed according to the seasons and weather conditions. For instance, particles collected in July, when the weather was calm, had δ values and C/N ratio representative of what is reported for plankton [START_REF] Peterson | Sulfur and carbon isotopes as tracers of saltmarsh organic matter flow[END_REF], while low δ 13 C values and high C/N ratio measured in June and October were characteristic at the studied site of surficial sediment resuspended after windy weather conditions [START_REF] Charles | Ecodynamics of PAHs at a peri-urban site of the French Mediterranean Sea[END_REF]. On these occasions, debris of vascular plants could be observed in the water samples and most photosynthetic pigments were degraded (data not shown). The sampling approach allowed thus to produce average values of δ and C/N that best reflected the food source for the shipworms and the co-occurring common filter-feeding species over the study period. It, however, should be noted that sampled SPOM particles ranged from 0.7 to 200 µm, and there is considerable heterogeneity in natural abundances of 13 C and 15 N among the discrete size classes of SPOM [START_REF] Rau | 13 C/ 12 C and 15 N/ 14 N variations among size-fractionated marine particles: implications for their origin and trophic relationships[END_REF]. δ 13 C and δ 15 N values of the smallest size fractions of SPOM can be as low as -25.3‰ and -0.5‰, respectively [START_REF] Rau | 13 C/ 12 C and 15 N/ 14 N variations among size-fractionated marine particles: implications for their origin and trophic relationships[END_REF]. At the study site, the contribution of picophytoplankton is maximal at periods of low total chlorophyll a biomass, that is to say in summer and fall, with contribution reaching more than 60% in September (Charles et al 2005). This, according to the findings of [START_REF] Rau | 13 C/ 12 C and 15 N/ 14 N variations among size-fractionated marine particles: implications for their origin and trophic relationships[END_REF], may to some extents explain the depletion in heavy isotopes measured in the SPOM fraction during the sampling, and so the gap between the bulk SPOM value and strict-filter feeding consumers which retain, sort and select particles differentially according to their size and nutritional quality [START_REF] Dubois | How picky can you be? Temporal variations in trophic niches of co-occurring suspension-feeding species[END_REF]Colombo, 2014, Cresson et al., 2016). As a consequence, even though this cannot dramatically change the general scheme of food source contributions, the comparison of particle sorting efficiency needs to be investigated to further understand interspecific variations in trophic niches among shipworm species [START_REF] Cresson | Feeding strategies of co-occurring suspension feeders in an oligotrophic environment[END_REF], Paalvast and van der Velde, 2013).

Conclusion

This study, involving alder wood, provides a first quantitative estimate of the contribution of food sources to the diet of B. carinata. Most of the assimilated carbon originated from the wood and N2 fixers contributed significantly to nitrogen requirements of the shipworm. However, it must be stressed that in a previous study conducted on Teredo navalis (Paalvast and van der Velde 2013), the overwhelming amount of carbon and nitrogen was derived from SPOM rather than from the wood source. Availability of food sources may change over time and space between studies, but stable isotope composition of wood also varies according to tree species and each species of teredinids has its own traits of life. As a consequence, the relative importance of the symbiotic association could then be better appreciated by testing, on the one hand, several tree and shipworm species, and, on the other hand, by performing 13

  N values, corresponding to values for invertebrate species (whole body) as reviewed by Caut et al. (2009).

  al 2009; Paalvast and van der Velde 2013, Yamaka et al., 2015). δ 15 N values of wood are more variable, showing particularly negative valuesas previously reported at this site by[START_REF] Carlier | Benthic trophic network in the Bay of Banyuls-sur-Mer (northwest Mediterranean, France): An assessment based on stable carbon and nitrogen isotopes analysis[END_REF] or for terrestrial plants[START_REF] Fogel | Unusually negative nitrogen isotopic compositions (δ 15 N) of mangroves and lichens in an oligotrophic, microbially-influenced ecosystem[END_REF]

Fig. 1 .

 1 Fig. 1. Carbon and nitrogen isotope compositions of potential food sources and consumers. The broken line indicates the δ 15 N average value expected for products resulting from the fixation of atmospheric nitrogen. The ellipse indicates the range of theoretical expected δ 15 N and δ13 C values for consumers of alder wood, considering trophic enrichments of 0.30‰ and 2.5‰ for δ13 C for δ 15 N, respectively. Dotted lines account for trophic shift of carbon and nitrogen.

Fig. 2 .

 2 Fig. 2. Bayesian mixing model dietary analyses for carbon (A) and nitrogen (B). Boxplots with 50 (in dark gray), 75 (in medium gray) and 95% (in light gray) credible intervals representing the proportion of Alnus glutinosa: Wood, suspended particulate organic matter: SPOM, and dissolved nitrogen: N2 in the diet of the shipworm Bankia carinata.

Figure S1 :

 S1 Figure S1: Experimental setup. Image of the alder wood packs at the study site within the Bay of Banyuls-sur-Mer, Mediterranean Sea, France.

Fig. 1 .

 1 Fig. 1. Carbon and nitrogen isotope compositions of potential food sources and consumers. The broken line indicates the δ 15 N average value expected for products resulting from the fixation of atmospheric nitrogen. The ellipse indicates the range of theoretical expected δ 15 N and δ 13 C values for consumers of alder wood, considering trophic enrichments of 0.30‰ and 2.5‰ for δ 13 C for δ 15 N, respectively. Dotted lines account for trophic shift of carbon and nitrogen.

Fig. 2 .

 2 Fig. 2. Bayesian mixing model dietary analyses for carbon (A) and nitrogen (B). Boxplots with 50 (in dark gray), 75 (in medium gray) and 95% (in light gray) credible intervals representing the proportion of Alnus glutinosa: Wood, suspended particulate organic matter: SPOM, and dissolved nitrogen: N2 in the diet of the shipworm Bankia carinata.

Table 1 .

 1 Cand15 N labelling experiments to quantitatively estimate the C and N pathways within the different potential components of the shipworms diet. Interfaces Littorales) topics, under the title: "De la forêt à la mer : Transfert et recyclage du bois mort sur la marge continentale méditerranéenne". bacteria enable a novel digestive strategy in a wood-feeding mollusk. Proceedings of the National Academy of Sciences, 111, E5096-E5104. https://doi.org/10.1073/pnas.1413110111 37. Paalvast, P., van der Velde, G., 2013. What is the main food source of the shipworm (Teredo navalis)? A stable isotope approach. Journal of Sea Research 80, 58-60. https://doi.org/10.1016/j.seares.2013.03.003 38. Parnell, A.C., Inger, R., Bearhop S., Jackson, A.L., 2010. Source partitioning using stable Stable isotope compositions and C/N ratios of potential food sources; n, number of samples. SPOM: Suspended particulate organic matter, SD: Standard deviation.

	isotopes:	coping	with	too	much	variation.	PLoS	ONE	5,	e9672.

Table 2 .

 2 Stable isotope compositions and C/N ratios of consumers. Feeding modes are: Xylo, xylotrophy, SF, suspension-feeding; n, number of individuals; SD: Standard deviation.

Table 1 .

 1 Stable isotope compositions and C/N ratios of potential food sources; n, number of samples. SPOM: Suspended particulate organic matter, SD: Standard deviation.

		Date or exposure
	Sources	duration in	n	 15 N (‰)	 13 C (‰)	C/N
		months
	SPOM	June 9	3	2.4 ± 0.1	-24.5 ± 0.1	6.5 ± 0.2
		July 11	4	4.8 ± 0.5	-22.6 ± 0.1	5.8 ± 0.6
		Sept 12	3	3.9 ± 0.5	-23.7 ± 0.1	6.2 ± 0.8
		Oct 24	4	2.6 ± 0.3	-24.3 ± 0.7	9.6 ± 2.7
		Time integrated		3.4 ± 1.1	-23.8 ± 0.8	7.0 ± 1.7
	Alder	0	5	-2.8 ± 0.7	-28.8 ± 0.5	158 ± 9
		2	4	-2.7 ± 0.1	-29.8 ± 0.1	195 ± 11
		4.5	4	-1.3 ± 0.5	-29.8 ± 0.3	157 ± 11
		Time integrated		-2.2 ± 0.8	-29.3 ± 0.5	170 ± 22

Table 2 .

 2 Stable isotope compositions and C/N ratios of consumers. Feeding modes are: Xylo, xylotrophy, SF, suspension-feeding; n, number of individuals; SD: Standard deviation.

	Class species	Feeding type	n	 15 N (‰)	 13 C (‰)	C/N
	Bivalvia
	Bankia carinata	Xylo/SF	25	0.5 ± 0.4	-27.5 ± 0.6	8.3 ± 2.5
	Mytilus galloprovincialis	SF	1	4.4 ± n.d.	-20.3	3.5 ± n.d.
	Hiatella arctica	SF	2	4.0 ± 0.07	-20.1 ± 0.2	3.8 ± 0.2
	Musculus subpictus	SF	1	4.1± n.d.	-20.8	4.0± n.d.
	Polychaeta
	Serpula vermicularis	SF	4	4.9 ± 0.1	-20.1 ± 0.3	3.9 ± 0.1
	Spirobranchus triqueter	SF	4	4.84± 0.2	-19.8 ± 0.2	3.9 ± 0.1
	Ascidiacea
	Ascidia conchilega	SF	1	4.6 ± n.d.	-21.3 ± n.d.	4.9 ± n.d.

Acknowledgments

We thank Gaël Guillou for stable isotope analysis. We are grateful to Joseph Garrigue, curator of the Réserve Naturelle Nationale de la Forêt de la Massane for providing alder wood for the field experiment.

We expressly thank Jean-Claude Roca and Bruno Hesse for the installation and recovery of wooden blocks at sea. We acknowledge Jan Ellis and Christian Nozais for the 'spell checking' and 'grammar checking' of the manuscript. This work was supported by the EC2CO (Écosphère Continentale et Côtière) national program of the CNRS-INSU, as part of the DRIL (Dynamique et Réactivité des


Stephane Brisset
email: stephane.brisset@centralelille.fr

Tuan-Vu Tran

PARETO-BASED BRANCH AND BOUND ALGORITHM FOR MULTIOBJECTIVE OPTIMIZATION OF A SAFETY TRANSFORMER

Keywords: Benchmark, Branch and bound, Multiobjective optimization, Pareto-optimality, Safety transformer

A multiobjective branch and bound method is presented and applied to the bi-objective combinatorial optimization of a safety transformer. New criteria are proposed for the branching and discarding. They are based on the Pareto dominance and contribution metric. The comparison with exhaustive enumeration and non-dominated sorting genetic algorithm confirms the solutions. It appears that exact and approximate methods are both very sensitive to their control parameters.

INTRODUCTION

The design of electromagnetic devices is mainly expressed in the literature in term of problem with continuous parameters. However, these problems are in the second part of the design process and often limited to the fine tuning of some parameters corresponding to the structure selected in the first part. Despite recent progress in topological [START_REF] Stolpe | Topology optimization -Recent applications and generalizations[END_REF] and combinatorial [START_REF] Amoiralis | A parallel mixed integer programming-finite element method technique for global design optimization of power transformer[END_REF]- [START_REF] Sourd | A multiobjective branch-and-bound framework: Application to the biobjective spanning tree problem[END_REF] optimizations, there is a lack of decision tools for the choice of the structure and materials when dealing with conflicting goals. At this stage, the parameters are mainly discrete and not sorted.

Moreover, the production in very small series practiced by some small and medium firms is supported by standards. It is thus a question of choosing among a great but finite number of solutions rather than to optimize some dimensions finely.

Optimization with discrete variables requires different concept than the conventional continuous one. The computation time of combinatorial optimization is also far more expensive. This is worsening in the design of electromagnetic devices because models are non-linear and time-consuming. Heuristic, Tabu search [START_REF] Glover | Modern Heuristic Techniques for combinatorial problems[END_REF], and Branch and Bound algorithm [START_REF] Amoiralis | A parallel mixed integer programming-finite element method technique for global design optimization of power transformer[END_REF]- [START_REF] Sourd | A multiobjective branch-and-bound framework: Application to the biobjective spanning tree problem[END_REF] can solve combinatorial problems. The formers compute approximate solutions in an affordable time while the latter find exact solutions with higher computing cost.

The first part of the paper is devoted to introduce the context of combinatorial optimization in electrical machines and the main issues for solving this kind of problems. In the second part, the mechanisms of BB algorithm are explained and new criteria for the branching and the initialization are proposed for multiobjective problems. The Pareto-based branch and bound algorithm is applied to the bi-objective optimization of a safety transformer. Results are compared to exhaustive enumeration and non-dominated sorting genetic algorithm. Finally, some conclusions and prospects are given.

CONTEXT AND ISSUES

Context

In the design of electromagnetic devices and especially electrical machines, the variables in optimization problems are often discrete. The variables can be integer and sorted like the number of slots and the number of magnets. When dealing with manufacturing in small series, continuous variables like the diameter of wire and the dimensions of magnetic core may be changed to discrete variables with values taken in the manufacturers' catalogs to take advantage of the price and availability.

Other design variables are more difficult to handle in optimization. They are discrete and cannot be sorted. For instance, the topologies of a rotating electrical machine are inner rotor, outer rotor, and axial flux as shown in Figure 1. The materials are also impossible to sort. This is illustrated by the case of conducting materials. If we compare copper to aluminum, copper has less losses but aluminum is lighter. So they cannot be sorted if both criteria appear in the optimization problem.

Issues

Three of the main issues that arise with combinatorial optimization are detailed to highlight the difficulties encountered in solving such problems. The first one is called integrality gap and is illustrated in the mathematical example with two design variables, one objective to minimize and two inequality constraints. The optimization problem is expressed in Eq. ( 1). One simple idea is to search the solution of combinatorial problem in the neighborhood of the solution of the relaxed problem expressed in Eq. ( 2) that is the problem with all integer variables changed to continuous. This last is the blue disc, written 𝑋 𝑐 * in Figure 2. Among the four integer solutions in its vicinity, the best one is the blue square 𝑋 𝑑 . It can be seen that the distance with the integer optimum 𝑋 𝑑 * at red square is significant. The second issue is that derivatives exist only with continuous variables. Therefore this information cannot be used to define a search direction to guide the optimization. As a consequence, the number of evaluations required for combinatorial optimization is very high.

min 𝑥 1 ,𝑥 2 𝑓(𝑥 1 , 𝑥 2 ) = (𝑥 1 -3) 2 + (𝑥 2 -3) 2 with 𝑥 1 ∈ {0.5, 1.0, 1.5, 2.0, 2.5} (
For multiobjective optimization, scalarization techniques cannot be applied to find the complete set of solutions. It is proved that the weighted sum of objectives fails to find it [START_REF] Ehrgott | Multicriteria Optimization[END_REF] and the authors' experience is that epsilon-constraint is inefficient and unable to find the complete set of non-dominated solutions.

BRANCH AND BOUND ALGORITHM

The principle of branch and bound algorithm is to divide the search domain into subdomains with lower dimensionality. The separation is made according to hyperplanes orthogonal to the variables axis. A tree spanning all combinations is erected (Figure 3). Given n, the number of integer design variables, there are n levels in the tree and one root at ground level. At the root node, all integer variables are relaxed, i.e. transformed to continuous variables. At the next levels, the evaluation of partial solutions at the nodes is an optimization with some integer variables relaxed and others being constant. Considering level k, the first k integer variables are fixed to constant values and the n-k lasts are relaxed. Therefore, the dimensionality of the optimization problem for the evaluation of partial solutions is decreasing while the level number increases. The leaves are at the last level with all integer variables set to constant value.

Figure 3. Spanning tree of combinations used in BB algorithm

The mechanisms used in BB algorithm are bounding, branching, and pruning as detailed below.

Branching scheme

Branching is the mechanism used to explore the tree in order to find all the solutions in a reduced computing time and memory requirement. The branching schemes available for BB are depth-first, best-first, and hybrid.

The former scheme is used here. This way, the tree is explored recursively by following the branch leading to the most promising node among all branches connected to the current node. The descent is continued until a leaf is reached. In order to reduce the amount of memory used, the tree is not described entirely during the initialization of algorithm but expanded progressively by adding a sub-tree after each branching. The depth-first scheme also offers the advantage to be the fastest to reach leaves and to enable the pruning of the tree.

The definition of the most promising node is based on a criterion to compare nodes to each other. In BB algorithms, this criterion is the bound.

Bounding method

BB methods have proved to perform well on many combinatorial optimization problems on the condition to provide a good bounding function. The bound can be computed by using an optimization process on a subproblem with some discrete variables fixed to feasible values and others relaxed, i.e. transformed to continuous variables. This makes sense for discrete variables that can be sorted if the objective and constraint functions could be defined for any real value within the range bounded by the minimum and maximum discrete values.

For the sake of simplicity, we consider first a single objective optimization as stated by the single axis in Figure 4. As the search sub-domain of nodes at the level k+1 is included in the domain at level k, the minimum of the problem for a node at level k+1 is higher or equal to the one for the node at level k to whom it is connected. Therefore, the partial solution at level k stands as a lower bound for all nodes at higher levels connected to it. Moreover, if no feasible solution is found at level k, no solution could be found at next levels. As a consequence, the most promising node is the partial solution with the lowest objective value and all constraints fulfilled. where 𝑣 𝑖 is the number of values allowed for the discrete variable 𝑥 𝑖 . The BB algorithm explores the tree iteratively by starting a descent each time that the leaf level is reached until no node remains. So, this is important to prune the tree as soon as possible to reduce the number of iterations and evaluations.

It is reminded that the partial solution is a lower bound. Therefore, a node is discarded if the value of the objective function is above the one of the best leaf as illustrated in Figure 4.

Example with two design variables

In Figure 5 is the spanning tree for the example with two variables expressed in Eq. ( 1). At level 1, the domain is separated according to 𝑥 1 . The five partial solutions are found by solving optimization problems with 𝑥 1 constant and 𝑥 2 continuous. During the first descent, i.e. iteration, the branching is made to the node number 4 with lowest objective value. At level 2, the subdomain corresponding to 𝑥 1 = 2.0 is separated according to 𝑥 2 . All nodes are leaves that are evaluated by a computation of objective and constraints. Four leaves are infeasible. The first descent is completed and the objective value of the best leaf is 7.25.

The second iteration starts with the elimination process. The node number 1 is discarded because its objective value (8.982) is higher than best leaf one. Three nodes remain and the descent starts from the node with the lowest objective value (6.014) whose number is 3. A new best leaf is found at leaf number 12 with an objective value of 6.25.

All remaining nodes have higher objective value than the new best leaf and are discarded. The algorithm stops because no better solution can be found. For this reason, BB algorithm is called exact method.

MULTIOBJECTIVE BRANCH AND BOUND ALGORITHM

The branching scheme is kept unchanged for multi-objective optimization branch and bound (MOBB) algorithm and new criteria are proposed to define the most promising node such as the elimination strategy.

Bounding method

Various approaches can be considered to deal with multiple conflicting goals. On one hand, it is tempting to reduce the number of objective to a single one to ease the comparison of partial solutions as in mono-objective BB algorithm. Several methods have been reported in the literature for MOBB algorithms such as scalarization (weighted sum of objectives, etc.).

On the other hand, the solutions of a multi-objective optimization problem are non-dominated solutions and referred as Pareto front when plotted in the objective space. Several metrics are proposed in the literature to assess the performances of multi-objective algorithms and to compare Pareto fronts [START_REF] Ehrgott | Multicriteria Optimization[END_REF]- [START_REF] Joshua | A tutorial on the performance assessment of stochastic multi-objective optimizers[END_REF]:

-D-metric assess the convergence to a reference Pareto front, --metric assess uniformity of points along the Pareto front, --metric compare Pareto front in term of extent, -Coverage metric evaluates the dominated area given by a front, and -Contribution metric evaluates the proportion of Pareto solutions given by each front. Intermediate approaches in terms of computing time and precision are reported in literature, for instance by using the ideal point whose coordinates are minimum values of objectives minimized separately [START_REF] Sourd | A multiobjective branch-and-bound framework: Application to the biobjective spanning tree problem[END_REF].

In this paper, partial solutions are compared by using the Pareto dominance and the bound is computed with the contribution metric [START_REF] Erbas | Multiobjective optimization and evolutionary algorithms for the application mapping problem in multiprocessor system-on-chip design[END_REF]. Therefore, the bound is no more an absolute value but a relative measurement of a node contribution among all the nodes at the same level linked to the same node at the previous level. This bounding method is valid for the depth-first branching scheme only. The most promising node is the one with the highest contribution.

As illustrated in Figure 6, among the solutions nodes at level k+1, nine are non-dominated. Four points belong to pink node and five to green node. As this last has the highest contribution, it is selected for branching.

Figure 6. Branching and discarding in multi-objective BB algorithm

Elimination strategy

The set of non-dominated leaves allows defining a dominated region as shown in blue in Figure 6. It is reminded that the partial solution is a lower bound that dominates all solutions at the next levels. Therefore, if all solutions at one node are within the dominated region then the node is discarded.

Initialization

An initial approximation of the Pareto front can lead to reduce the number of evaluations required by the multi-objective branch and bound algorithm. This front may be found by metaheuristic methods but here, it is proposed to use the partial front at root node.

In the left part of Figure 7, the solutions of partial front at root node are plotted in the objectives space and in the right part they are plotted in the design variables space. For each solution, the integer solutions in its vicinity are evaluated and kept if all constraints are fulfilled. The set of non-dominated initial leaves becomes the initial front.

Mixed variables

BB algorithms can solve mixed variables problems, i.e. problems with discrete, integer, and continuous variables. Continuous variable can take any real value within a range. Integer variables are sorted and assumed to be relaxable. This means that the objectives and constraints functions are defined for any real value between the minimum and the maximum discrete values and are smooth. This can be critical when using a finite element model and considering an integer variable like the number of slots. Discrete variables cannot be sorted and are not relaxed during the partial evaluations of BB algorithm.

Therefore, a first spanning tree is raised with the discrete variables and its leaves become the roots of subtrees with integer and continuous variables. The leaves of the sub-trees are evaluated by an optimization where all the integer variables are set to constant value and the continuous variables are unchanged.

COMBINATORIAL OPTIMIZATION OF SAFETY TRANSFORMER

Device and models

The electromagnetic device studied here is a safety isolating transformer [START_REF] Tran | Proposal of a benchmark for multi-level optimization with 3D finite element model[END_REF] (Fig. 8). An analytical model is used to compute the objectives and constraints. This model accepts continuous parameters as inputs and can thus be used with relaxation techniques.

Figure 8. Safety transformer variables

Thermal and magnetic phenomena are both modeled by using 3D FEA on one eighth of the transformer due to the symmetries (Fig. 9). All magnetic and electric quantities are assumed sinusoidal. Full-load and no-load simulations are used to compute all the characteristics. The iron loss is computed with Steinmetz formula and the leakage inductances are calculated with the magnetic co-energy. The magneto-thermal coupling is weak and consistency loop requires 2 hours.

To reduce the computing time, a lumped-mass model is preferred. The additional modeling hypotheses are uniform distribution of induction in the iron core, no voltage-drop due to the magnetizing current, and uniform temperatures in coils and lamination. The multiphysic coupling is strong and the computing time is 50 ms.

Optimization problem

The problem contains 7 integer design variables that are shown in Fig. 8. There are three parameters (a, b, c) for the shape of the lamination, one for the frame (d), two for the section of enameled wires (S1, S2), and one for the number of primary turn n1. As the first four integer variables are linked, the problem becomes a mixedinteger problem with one discrete variable for the combination the firsts four and the three remaining integer variables. There are 24 types of lamination from catalogue r.bourgois®, 62 possible matches between the laminations EI and the frames from catalogue isolectra-martin®, and 62 types of enameled wires from invex®. The number of primary turn n1 is integer but only 1001 values are allowed, leading to 238,566,328 combinations.

There are 7 inequality constraints in this problem. The copper and iron temperatures Tcond, Tiron respectively are less than 120°C and 100°C. The magnetizing current I10/I1 and drop voltage V2/V2 are less than 10%, the filling factor of both coils f1, f2 is lower than 1, and the residue of coupled equations is less than 10 -6 . The objective is to minimize the total mass Mtot of iron and copper materials and to maximize the efficiency η:

           V V residue f f I I C T C T M iron cond tot    EI c d b a  , , , ( 62 configurations) (4)
  1200 , , 200 1   n (1001 values) W S W S   2 1
(62 values each)

Results

A time consuming exhaustive enumeration of solutions is performed with distributed computing on 24 cores and has found 749 non-dominated solutions with 238,566,328 evaluations of the model.

For MOBB, the bound is computed at each node with the Pareto set of the sub-problem. This set is found with sequential quadratic programming (SQP), -constraint, and multi-start. Its precision is depending on the number of starting points and the number of Pareto points requested. With 20 Pareto points and one starting point, the whole set of solutions is found. Initialization of MOBB with integer solutions in the vicinity of relaxed ones slightly decreases the number of evaluations by 7% as detailed in Table 1.

The number of evaluations and the size of the Pareto set found by MOBB are compared with NSGA-II [START_REF] Deb | A fast and elitist multiobjective genetic algorithm NSGA-II[END_REF]. An implementation of this algorithm is given by Song Lin in Matlab Central (NGPM v1.4, 2011). For NSGA-II, the results depend on the population size and the number of generations. One thousand individuals and two hundred generations are sufficient to find the whole solutions with five times less evaluations than MOBB. On Figure 10 are two Pareto fronts. The Pareto front found by solving the relaxed problem at root node is in blue. It is continuous and contains one hundred points linked by lines. The solutions of the combinatorial problem are in red and the front is discontinuous. All the solutions of the combinatorial problem are dominated by those of the relaxed one, as expected.

Discussions

The main weakness of MOBB is the computing burden to evaluate partial solutions while the main drawbacks of genetic algorithms that are the precision of the solution and especially the respect of constraints vanish in integer optimization. It could differ with a mixed-integer optimization problem because of the lack of precision of genetic algorithm for the continuous components of the solutions and the fulfillment of constraints when solutions are on the limit-state.

One track to reduce the number of evaluations of MOBB is to lower the precision of bounds as their purpose is to compare the partial solutions among them and not to calculate the integer solutions. However, we shall remain cautious because Figure 10 shows that the partial front at root node is close to the integer solutions at leaves while they are situated at the opposite extrema of the spanning tree. It may also be possible to be less sensitive to the number of Pareto points by using a most accurate criterion than contribution metric to compare partial solutions.

OSYCZKA AND KUNDU TEST FUNCTION

Optimization problem

The problem is given in Eq. ( 5) and ( 6) where all the variables are continuous within bounds [START_REF] Osyczka | A new method to solve generalized multicriteria optimization problems using the simple genetic algorithm[END_REF]. It is changed to a combinatorial optimization problem by using a step of 0.25 for all variables. Obviously, one objective and one constraint are non-convex. As a consequence, SQP, -constraint, and multi-start may fail to find exact partial solutions of sub-problems.

                                     x x x x x x x x x x x x x x x x x x f x x x x x f

Results

An exhaustive enumeration of solutions has found 73 non-dominated solutions with 497,954,225 evaluations. In table 2, it can be seen that MOBB fails to find the complete set of solutions in an affordable number of evaluations while NSGA2 succeed. This test highlights that the proposed algorithm has poor performances with non-convex functions. CONCLUSIONS MOBB and NSGA-II succeed to find the complete Pareto set of the bi-objective safety transformer combinatorial problem with integer variables but NSGA-II is yet four to five times faster. The performances of both algorithms are sensitive to their control parameters.

Figure 1 .

 1 Figure 1. Topologies of a rotating electrical machine

Figure 4 .

 4 Figure 4. Branching and discarding in single-objective BB algorithm Elimination strategy One of the advantages of the branching scheme called depth-first if to be the fastest to reach the leaves. This is done by using ∑ 𝑣 𝑖 𝑛-1 𝑖=1partial evaluations, i.e. optimizations and 𝑣 𝑛 evaluations while the number of combinations, i.e. leaves is ∏ 𝑣 𝑖 𝑛 𝑖=1

Figure 5 .

 5 Figure 5. Spanning tree for the example with two design variables

Figure 7 .

 7 Figure 7. Initial front for multi-objective BB algorithm

Figure 9 .

 9 Figure 9. Thermal and electromagnetic models of the safety transformer

Figure 10 .

 10 Figure 10. Pareto front found by of MOBB (discrete) and initial front (relaxed)

  On Figure11are two Pareto fronts. The Pareto front found by solving the relaxed problem is in blue. It is found by SQP and -constraint with 500 values of  and 500 starting points. The solutions of the combinatorial problem are in red.

Figure 11 .

 11 Figure 11. Pareto fronts for discrete and relaxed variables

𝑥 1 , 𝑥 2 ) ≤ 0 Figure 2. Mathematical

	1)	𝑥 2 ∈ {0.5, 1.0, 1.5, 2.0, 2.5}
	s.t.	𝑔 1 (𝑥 1 , 𝑥 2 ) = 10 𝑥 1 + 9 𝑥 2 3 -17 ≤ 0
		𝑔 2 (𝑥 1 , 𝑥 2 ) = 2 𝑥 1 + 5 𝑥 2 -8.3 ≤ 0
	min 𝑥 1 ,𝑥 2	𝑓(𝑥 1 , 𝑥 2 )
	with	0.5 ≤ 𝑥 1 ≤ 2.5
	(2)	0.5 ≤ 𝑥 2 ≤ 2.5
	s.t.	𝑔 1 (𝑥 1 , 𝑥 2 ) ≤ 0
		𝑔 2 (

example with two design variables (x1 is abscissa and x2 is ordinate)

Table 1 . Comparison of algorithms

 1

	algorithm	population size or	generations or	number of	size of Pareto set
		Pareto points	starting points	evaluations
	MOBB	10	1	592,784	732

Table 2 . Comparison of algorithms

 2

	algorithm	population size or	generations or	number of	size of Pareto set
		Pareto points	starting points	evaluations
	MOBB	100	10	37,576,081	47
	MOBB + initialization	100	10	17,356,443	47
	MOBB	200	10	85,519,229	60
	MOBB + initialization	200	10	35,677,027	57
	NSGA2	76	100	7,600	73

Three prospects are proposed to reduce the number of evaluations required by exact methods. The first one is to do a preliminary pruning of the spanning tree by evaluating explicit constraints on the filling factor of both coils.

The main drawback of the algorithm proposed is the high number of evaluations required to compute the partial solutions at nodes. The second prospect is to lower the precision of bounds in order to reduce the computation burden and to use a new criterion to compare partial solutions.

Branch and cut algorithms process is similar to BB but the separation is made according to cut plane defined by the constraints of the optimization problem. The third prospect is to test it to see if it is faster than BB for hard constrained problems.

Osyczka and Kundu test function highlights that the proposed algorithm has poor performances with nonconvex objective and constraint.