Aurélie Chauffour

Nacer Lounis

Koen Andries

Vincent Jarlier

Nicolas Veziris

Alexandra Aubry
email: alexandra.aubry@sorbonne-universite.fr

Minimal effective dose of bedaquiline administered orally and activity of a long acting formulation of bedaquiline in the murine model of leprosy Bedaquiline for the treatment of leprosy

Keywords: leprosy, bedaquiline, footpad mouse model

Background

Bedaquiline (BDQ), by targeting the electron transport chain and having a long half-life, is a good candidate to simplify leprosy treatment. Our objectives were to (i) determine the minimal effective dose (MED) of BDQ administered orally, (ii) evaluate the benefit of combining two inhibitors of the respiratory chain, BDQ administered orally and clofazimine (CFZ)) and (iii) evaluate the benefit of an intramuscular injectable long-acting formulation of BDQ (intramuscular BDQ, BDQ-LA IM), in a murine model of leprosy.

Methodology/ Principal Findings

To determine the MED of BDQ administered orally and the benefit of adding CFZ, 100 fourweek-old female nude mice were inoculated in the footpads with 5x10 3 bacilli of M. leprae strain THAI53. Mice were randomly allocated into: 1 untreated group, 5 groups treated with BDQ administered orally (0.10 to 25 mg/kg), 3 groups treated with CFZ 20 mg/kg alone or combined with BDQ administered orally 0.10 or 0.33 mg/kg, and 1 group treated with rifampicin (RIF) 10 mg/kg. Mice were treated 5 days a week during 24 weeks.

To evaluate the benefit of the BDQ-LA IM, 340 four-week-old female swiss mice were inoculated in the footpads with 5x10 3 to 5x10 1 bacilli (or 5x10 0 for the untreated control group) of M. leprae strain THAI53. Mice were randomly allocated into the following 11 groups treated with a single dose (SD) or 3 doses (3D) 24h after the inoculation: 1 untreated group, 2 treated with RIF 10 mg/kg SD or 3D, 8 treated with BDQ administered orally or BDQ-LA IM 2 or 20 mg/kg, SD or 3D.

Twelve months later, mice were sacrificed and M. leprae bacilli enumerated in the footpad.

orally against M. leprae in this model is therefore 3.3 mg/kg. The combination of CFZ and BDQ 10-fold lower than this MED did not significantly increase the bactericidal activity of CFZ. The BDQ-LA IM displayed similar or lower bactericidal activity than the BDQ administered orally.

Conclusion

We demonstrated that the MED of BDQ administered orally against M. leprae was 3.3 mg/kg in mice and BDQ did not add significantly to the efficacy of CFZ at the doses tested. BDQ-LA IM was similar or less active than BDQ administered orally at equivalent dosing and frequency but should be tested at higher dosing in order to reach equivalent exposure in further experiments.

Introduction

Leprosy remains a major health problem worldwide despite being one of the oldest infectious diseases, reported for more than 2000 years. The leprosy elimination goal as a public health problem set by the World Health Organization, aiming for a global prevalence rate of < 1 patient in a population of 10,000, was achieved in 2000, but up to 200,000 new cases are still reported each year [START_REF] Organization | Global leprosy strategy 2016-2020: accelerating towards a leprosy-free world -2016 operational manual[END_REF]. The worldwide use of leprosy drugs starting in the 1980s and their access at no cost for patients since 1995 were tremendous in the ability to achieve leprosy elimination [START_REF]Chemotherapy of leprosy for control programmes : report of a WHO study group[END_REF]. Nowadays, the WHO global strategy targets zero leprosy by 2030, but this goal has been hindered by sharp reduction of leprosy case detection during 2020-2021 due to the Covid-19 pandemic [START_REF][END_REF]. As with other bacteria of medical interest, antimicrobial resistance is observed in the causative agent Mycobacterium leprae in several parts of the world, despite multidrug therapy being the recommended standard leprosy treatment to avoid resistance selection since 1982. The first treatment of leprosy, consisting of a monotherapy of dapsone, led to the emergence of drug-resistance [START_REF] Pearson | Primary dapsone-resistant leprosy[END_REF]. Despite the addition of rifampicin (RIF) in the 1960s, drug-resistant strains quickly emerged [START_REF] Jacobson Robertr | RIFAMPIN-RESISTANT LEPROSY[END_REF].

Moreover, the length of the treatment leads to a poor compliance by patients and may favor the emergence of resistant strains. Therefore, to simplify and to facilitate the direct observation of treatment, a shorter, fully supervisable, monthly-administered multidrug regimen for leprosy is highly desirable [START_REF] Ji | Bactericidal activities of R207910 and other newer antimicrobial agents against Mycobacterium leprae in mice[END_REF]. Finally, in addition to patients whose M. leprae isolates are resistant to RIF, special regimens are also required for individual patients who cannot take RIF because of allergy, concomitant drug interaction or intercurrent disease such as chronic hepatitis.

In 2005, a newly discovered class of antibiotics, the diarylquinoline, was reported to be highly bactericidal against M. tuberculosis in mice and later in the mouse models for M. leprae [START_REF] Andries | A diarylquinoline drug active on the ATP synthase of Mycobacterium tuberculosis[END_REF][START_REF] Gelber | The Diarylquinoline R207910 Is Bactericidal against Mycobacterium leprae in Mice at Low Dose and Administered Intermittently[END_REF]. The lead compound, bedaquiline (BDQ), also called R207910 or TMC207, inhibits an enzyme belonging to the electron transport chain, the ATP synthase, by binding to the subunit c of the enzyme, leading to a decrease in bacterial metabolism. The bactericidal activity of BDQ administered orally against M. leprae observed in mice is similar to that of moxifloxacin and RIF supporting the launch of a clinical trial aiming at evaluating BDQ efficacy in multi bacillary (MB) leprosy [START_REF] Research | An Open-Label Study to Evaluate the Efficacy and Safety of TMC207 in Subjects With Multibacillary Leprosy[END_REF]. BDQ administered orally is currently the only new drug under clinical trial for leprosy treatment [START_REF] Cole | Massive gene decay in the leprosy bacillus[END_REF].

Interestingly, M. leprae does not possess all the proteins along the electron transport chain [START_REF] Cole | Massive gene decay in the leprosy bacillus[END_REF] suggesting that associating inhibitors acting at different enzymes belonging to it may act synergistically and may display strong bactericidal anti-leprosy activity. Clofazimine (CFZ), whose mechanism of action is not fully understood, targets the electron transport chain at the level of the menaquinone. CFZ has been shown to be effective against leprosy [START_REF] Yawalkar | Lamprene (clofazimine) in leprosy. Basic information[END_REF][START_REF] Black | Energy Metabolism and Drug Efflux in Mycobacterium tuberculosis[END_REF][START_REF] Ji | Bactericidal activities of combinations of new drugs against Mycobacterium leprae in nude mice[END_REF][START_REF] Ji | Clinical trial of ofloxacin alone and in combination with dapsone plus clofazimine for treatment of lepromatous leprosy[END_REF][START_REF] Gelber | The killing of Mycobacterium leprae in mice by various dietary concentrations of clofazimine and ethionamide[END_REF][START_REF] Levy | Activity of four clofazimine analogues against Mycobacterium leprae[END_REF][START_REF] Holmes | Effect of rifampin, clofazimine, and B1912 on the viability of Mycobacterium leprae in established mouse footpad infection[END_REF][START_REF] Shepard | Combinations involving dapsone, rifampin, clofazimine, and ethionamide in the treatment of M. leprae infections in mice[END_REF][START_REF] Shepard | Minimal effective dosages in mice of clofazimine (B,663) and of ethionamide against Mycobacterium leprae[END_REF].

The anti-leprosy activity of the association of the two drugs acting on electron transport chain, e.g. BDQ and CFZ, deserve to be evaluated. In addition, the genes encoding the M. tuberculosis MmpS5-MmpL5 efflux pump and repressor (mmpR, Rv0678) whose mutations are implicated in BDQ resistance are absent in M. leprae where it might contributes to a higher potency of bedaquiline against M. leprae compared to M. tuberculosis [START_REF] Cole | Massive gene decay in the leprosy bacillus[END_REF][START_REF] Melly | MmpL Proteins in Physiology and Pathogenesis of M. tuberculosis[END_REF].

Two key properties of drugs administered in LAI (Long Acting Injectable) formulations are low aqueous solubility to preclude the rapid dissolution and release of the active drug substance, and a reasonably long pharmacokinetic (PK) elimination half-life, i.e., slow clearance from the body. For an antimicrobial, another desired property is high potency, negating the need for high concentrations in the blood and allowing lower drug doses to be injected. Bedaquiline seems to be highly potent against M. leprae [9] at a lower dose than against M. tuberculosis. It has high lipophilicity (logP, 7.3), and a long half-life (about 24 h, functionally or effectively) which makes it suitable for use in an LAI formulation [START_REF] Kaushik | Activity of a Long-Acting Injectable Bedaquiline Formulation in a Paucibacillary Mouse Model of Latent Tuberculosis Infection[END_REF]. The efficacy of LAI BDQ has been already demonstrated in a latent tuberculosis infection mouse model [START_REF] Kaushik | Activity of a Long-Acting Injectable Bedaquiline Formulation in a Paucibacillary Mouse Model of Latent Tuberculosis Infection[END_REF][START_REF] Kaushik | Efficacy of Long-Acting Bedaquiline Regimens in a Mouse Model of Tuberculosis Preventive Therapy[END_REF][START_REF] Nguyen | Linking In Vitro Intrinsic Dissolution Rate and Thermodynamic Solubility with Pharmacokinetic Profiles of Bedaquiline Long-Acting Aqueous Microsuspensions in Rats[END_REF]. Due to the very slow doubling time of M. leprae, a unique administration of LAI BDQ could also be considered.

In our work, we aimed to (i) determine the minimum effective dose (MED) of the BDQ administered orally, (ii) evaluate the benefit of combining BDQ administered orally and CFZ, and (iii) evaluate the benefit of a BDQ long-acting formulation in a murine model of leprosy.

Methods

Ethics statement

The experimental project was favorably evaluated by the ethics committee n°005 Charles Darwin localized at the Pitié-Salpêtrière Hospital. Clearance was given by the French Ministry

Materials

In both experiments, mice were infected with a M. leprae THAI53 strain. This strain was fully susceptible to the common antileprosy drugs (i.e. RIF, dapsone, CFZ and fluoroquinolones) [START_REF] Matsuoka | The history of Mycobacterium leprae Thai-53 strain[END_REF]. The suspension used to inoculate mice was prepared from mice already infected by this isolate one year earlier. Shepard and Mac Rae's method was used to prepare the suspension [START_REF] Shepard | A method for counting acid-fast bacteria[END_REF]. Briefly, the tissue from the footpads was aseptically removed and then grinded by using a GentleMacs Octo Dissociator (Miltenyi) under a volume of 2 ml of Hanks' balanced salt solution. Ten µl of the prepared suspension were taken to perform a Ziehl-Neelsen staining to count M. leprae Acid Fast Bacilli (AFB). Suspensions needed to inoculate mice were then further diluted in Hanks' balanced salt solution.

Respectively four-week-old nude (NMRI-Foxn1 nu/nu ) for the determination of the minimal effective dose of BDQ administered orally and swiss mice for the evaluation of the BDQ-LA IM were purchased from Janvier Labs, Le Genest Saint Isle, France. The nude mice model aims to mimic low-immunity leprosy and the swiss, high-immunity leprosy; therefore those two models will respectively represent multibacillary vs paucibacillary leprosy. RIF and CFZ were purchased from Merck, France; BDQ administered orally and BDQ-LA IM were provided by Johnson and Johnson, Belgium.

Infection of mice with M. leprae and treatment

First experiment: determining the minimal effective dose of BDQ administered orally against M. leprae and the contribution of CFZ when combined with BDQ. We adapted the continuous method to determine the MED of BDQ administered orally [START_REF] Shepard | Effect of Several Anti-Leprosy Drugs on Multiplication of Human Leprosy Bacilli in Foot-Pads of Mice[END_REF]. The MED is defined as the lowest dose of a drug that inhibits the growth of M. leprae, i.e. corresponding to the group where all the mice footpads remain negative. One hundred 4-week-old female nude mice were infected in the left hind footpad with 0.03 ml of the M. leprae isolate THAI53 according to Shepard's method [START_REF] Shepard | The experimental disease that follows the injection of human leprosy bacilli into foot-pads of mice[END_REF] with an inoculum of 5x10 3 AFB/ footpad. Mice were then randomly allocated into one untreated control group and 9 treated groups of 10 mice each: RIF 10 mg/kg, BDQ administered orally 0.10, 0.33, 1, 3.3 or 25 mg/kg, CFZ 20 mg/kg and combinations of BDQ administered orally 0.10 or 0.33 mg/kg and CFZ 20 mg/kg.

Treatment was given one month after inoculation, five days a week during 24 weeks by oral gavage under a volume of 0.2 ml per mouse. The RIF and CFZ dosages used in our study are the current doses used in the murine model of leprosy [START_REF] Ji | Bactericidal activities of combinations of new drugs against Mycobacterium leprae in nude mice[END_REF][START_REF] Gelber | Detection of persisting Mycobacterium leprae by inoculation of the neonatally thymectomized rat[END_REF].

Second experiment: comparing the bactericidal activity of BDQ administered orally and BDQ-LA IM against M. leprae. We used the proportional bactericidal method that allows to measure the bactericidal activity of a compound [START_REF] Colston | The "proportional bactericidal test": a method for assessing bactericidal activity in drugs against Mycobacterium leprae in mice[END_REF]. Three hundred and forty 4-week-old female swiss mice were infected in the left hind footpad with 0.03 ml of the M. leprae isolate according to Shepard's method [START_REF] Shepard | The experimental disease that follows the injection of human leprosy bacilli into foot-pads of mice[END_REF]. Mice were inoculated with three different inocula of 5x10 3 , 5x10 2 , 5x10 1 AFB/ footpad except for the untreated control which was also inoculated with one 5x10 0 additional group. Mice were randomly allocated into one untreated control group and 10 treated groups of 10 mice each: RIF 10 mg/kg, BDQ administered orally 2 or 20 mg/kg, BDQ-LA IM 2 or 20 mg/kg. Treatment was given by oral gavage under a volume of 0.2 ml per mouse, except for the BDQ-LA IM which was injected intramuscularly under a volume of 0.012 ml per thigh and both thighs were injected at the same time. Treatment for all drugs was given as a single (SD) or three doses (3D) for BDQ administered orally 2 or 20 mg/kg and BDQ-LA IM 2 or 20 mg/kg and began the day after inoculation for the SD, and 4 and 8 weeks later for 3D.

Assessment of the effectiveness of the treatment

To permit multiplication of M. leprae to a detectable level, mice were held 12 months in the animal facility. Mice were then euthanized and tissues from their footpad were removed aseptically and homogenized under a volume of 2 ml of Hank's balanced salt according to the Shepard's method [START_REF] Shepard | The experimental disease that follows the injection of human leprosy bacilli into foot-pads of mice[END_REF]. M. leprae bacilli were considered to have multiplied (i.e. survived the treatment) if those footpads were found to contain 10 5 acid-fast bacilli, regardless of the size of the inoculum.

Statistical analysis

First experiment: determine the MED of BDQ administered orally against M. leprae. A Fisher exact test was performed. A p-value <0.05 was considered to be statistically significant by standard evaluation. For multiple comparisons between the groups, Bonferroni's correction was applied, i.e., the difference would be significant at the 0.05 level only if the P value adjusted to the number of groups: 0.05/n in which n was defined as the number of primary comparisons. Thus, the corrected P was 0.05/10 = 0.005.

Second experiment: compare the bactericidal activity of BDQ administered orally and BDQ-LA IM against M. leprae. The proportion of viable M. leprae after treatment was determined from the infectious dose required to show multiplication in 50% of the inoculated mice. The significance of the differences between the groups was calculated by the Spearman and Kärber method [START_REF] Shepard | Statistical analysis of results obtained by two methods for testing drug activity against Mycobacterium leprae[END_REF]. A p-value <0.05 was considered statistically significant by standard evaluation. For multiple comparisons between the groups, Bonferroni's correction was applied, i.e., the difference would be significant at the 0.05 level only if the P value adjusted to the number of groups: 0.05/n in which n was defined as the number of primary comparisons. Thus, the corrected P was 0.05/11 = 0.0045.

Results

Minimal effective dose of BDQ administered orally (Table 1, Figure 1)

After one year of observation, all footpads were positive in the untreated control group (mean of 8.13±0.26 log 10 AFB per footpad), confirming the multiplication of M. leprae.

As compared to untreated control, RIF reduced the bacillary load since 4 mice had negative footpads after treatment (p=0.08).

All mice footpads remained positive after treatment in the 3 lowest BDQ administered orally doses (0.10, 0.33 and 1 mg/kg) but they reduced the bacillary load by 1 log10 AFB as compared to the untreated control (p=0.002, p=0.00001, p=0.0009 respectively). These 3 BDQ administered orally doses were also less bactericidal than RIF (10 mice remained positive after treatment in the 3 BDQ administered orally doses groups vs 6 mice positive in the RIF group) but it was not statistically significant (p=0.08). On the other hand, all the footpads were negative in the groups treated with the 2 highest BDQ administered orally doses 3.3 mg/kg and the 25 mg/kg which was statistically significant regarding untreated control (p=0.00001 for both groups). These groups displayed higher bactericidal activity than RIF 10 mg/kg but the difference was not statistically significant with Bonferroni's correction (p=0.01 for both groups). CFZ 20 mg/kg was bactericidal as compared to untreated control (p=0.03) even if 3 mice remained positive after treatment. CFZ was as bactericidal as RIF with a p value not statistically different (p=0.370). The combination of the 2 lowest BDQ administered orally doses (i.e. 0.10 or 0.33 mg/kg) with CFZ lead to a decrease in the number of AFB positive footpads as compared with CFZ but the reduction was not statistically significant (p=1 and 0.58, respectively).

Comparison of BDQ administered orally and BDQ-LA IM (Table 2)

From results in the untreated group, the proportion of viable M. leprae was estimated to be 2.75% of the total number of AFB inoculated. The percentage of viable bacilli killed under treatment ranged between 74.86% and 99.85% depending on the treatment group.

Compared to untreated control, the percentage of viable bacilli was smaller in all groups except BDQ 2 mg/kg SD either oral or LA IM; but the difference between the untreated group and the RIF 10 mg/kg SD was not significant after Bonferroni's adjustment. RIF 3D was significantly different from RIF SD (p<0.001).

When comparing groups treated with 2 mg/kg BDQ administered orally or BDQ-LA IM at the same frequency (SD or 3D), a similar number of footpads remained positive in both groups (p>0.05). When comparing groups treated with 20 mg/kg BDQ administered orally or BDQ-LA IM at the same frequency (SD or 3D), more footpads remained positive in the groups treated with BDQ-LA IM and the difference were statistically significant but not after Bonferroni adjustments whatever the number of administrations between BDQ and BDQ-LA IM groups (p>0.0045) (Table 2). a treatment was given by oral gavage, except for BDQ-LA IM which was administered intramuscularly b M. leprae bacilli were considered to have multiplied if the harvest from a footpad yielded >10 5 acid-fact bacilli c the proportion of viable M. leprae surviving the treatment could be calculated by estimating the "most probable number" of viable organisms. However, the estimation of the MPN is based on the assumption that the organisms are distributed randomly in an inoculum; in the case of M. leprae, this assumption is probably untenable, therefore, the preferred alternative is to calculate the "median infectious dose (ID50)", i.e. the number of organisms required to infect 50% of the mice as allowed by the Spearman-Kärber method (it requires that the titration be carried out over a range of 100% to 0%). In immunocompetent mice, if the largest inoculum is 5X10 3 M. leprae per footpad, a proportion of viable M. leprae as small as 0.00006 may be measured, then it is possible to calculate the proportion of viable M. leprae killed by the treatment by comparing the proportions of viable in tested and control mice. The significance of the differences between the groups was calculated by the Spearman and Kärber method [START_REF] Shepard | Statistical analysis of results obtained by two methods for testing drug activity against Mycobacterium leprae[END_REF] d calculated from the comparison of the proportion of viable organisms between untreated controls and the treated group.

e each treated group was compared to the untreated group of mice f drug was given under a single dose the day after infection g drug was given three times (the day after infection, 4 and 8 weeks later) h p value corresponds to the comparison of the 2 groups at the beginning of the line: RIF SD vs 3D, BDQ oral vs LA-IM at equivalent dosing and frequency fold lower MED of BDQ administered orally against M. leprae compared to that of M.

tuberculosis.

It's important to mention that a single dose of 2 mg/kg BDQ administered orally (total dose = 2 mg/kg) in experiment 2 was able to kill 75% of the viable M. leprae bacilli while 6 months of 5 days per week treatment of 1 mg/kg (total dose = 120 mg/kg) in experiment 1 was not able to kill any bacilli present at start of treatment but was only able to slow down the replication of the bacilli when compared to untreated mice. The main differences between experiments 1 & 2 are the use of immunocompromised mice and initiation of treatment one-month post-infection in experiment 1 and the use of immunocompetent mice and initiation of treatment the day after infection in experiment 2. The use of immunocompromised mice without T cell immunity allows a much better replication of M.

leprae and a much higher rate of viable bacilli while the replication rate and the rate of viable bacilli are much lower in immunocompetent mice which partially explain the differences in the effectiveness of the 2 regimens.

The treatment of leprosy needs to be based on a drug-combination to avoid the emergence of resistant-strains. BDQ targets the electron transport chain, which is also the target of the classical anti-leprosy drug CFZ [START_REF] Black | Energy Metabolism and Drug Efflux in Mycobacterium tuberculosis[END_REF][START_REF] Ji | Bactericidal activities of combinations of new drugs against Mycobacterium leprae in nude mice[END_REF]. Combining low doses of BDQ that were ineffective alone with CFZ enabled to reduce the numbers of AFB positive footpads (Table 1), however, BDQ did not add significantly to the efficacy of CFZ at the tested doses . These 2 combinations were slightly more bactericidal than RIF. Despite being not statistically significant, this result may suggest that a combination of drugs targeting the electron transport chain may be highly bactericidal against leprosy. Further experiments should be designed and performed in nude mice rather than in Swiss, since the much larger numbers of M. leprae viable organisms in nude mice than in immunocompetent Swiss mice may permit more accurate differentiation among various levels of bactericidal activities.

One of the main characteristics of M. leprae is its long doubling time (i.e. 14 days) suggesting an active drug with a long half-life would be a good choice against this bacterium. The drug must be effective at a low dose and under an intermittent administration, conditions that are currently sought for the treatment of leprosy. A new formulation of BDQ, called BDQ-LA IM, was found to be active in a mouse model of latent tuberculosis [START_REF] Kaushik | Efficacy of Long-Acting Bedaquiline Regimens in a Mouse Model of Tuberculosis Preventive Therapy[END_REF]. Despite the potential challenges of introducing an injectable formulation in the field, its improved pharmacokinetics properties compared to BDQ administered orally may allow reduction in the duration of the treatment and therefore increase adherence of patients to the treatment. The BDQ-LA IM formulation was tested at two dosages, 2 or 20 mg/kg, and with two frequencies (one, SD or three doses, 3D). Our results showed that at equivalent dosing and frequency, BDQ administered was similar, or even more active than LA-IM in our murine model of leprosy (Table 2). A possible explanation why the 20 mg/kg was less effective when administered as a LA injection rather than orally is that the doses tested in our experiment were too low for LA-IM BDQ. Indeed, in the studies evaluating the BD-LA IM in tuberculosis higher doses were used and it was shown that BDQ LA-IM 160 mg/kg generated a C max equivalent to that of 30 mg/kg BDQ administered orally. In support of this hypothesis is the fact that in the present study when comparing 20 mg/kg LA-IM BDQ to 2 mg/kg BDQ administered orally at equivalent dosing, the injectable form was more active than the oral form. Consequently BDQ LA-IM should be further tested at higher dosing against leprosy.

In conclusion, we found that the MED of BDQ administered orally against M. leprae was 3.3 mg/kg and that the combination of CFZ and BDQ may improve the bactericidal activity. BDQ-LA IM was similar or less active than BDQ administered orally at equivalent dosing and frequency but should be tested at higher dosing in order to reach equivalent exposure in further experiments. These findings open the path to the design of shorter BDQ-based treatment deserving to be evaluated in leprosy.

of

  Higher Education and Research under the number APAFIS#9575-2017030114543467 v3. Our animal facility received the authorization to carry out animal experiments (license number D75-13-08). The persons who carried out the animal experiments had followed a specific training recognized by the French Ministry of Higher Education and Research and follow the European and the French recommendations on the continuous training.

Figure 1 .

 1 Figure 1. Multiplication of M. leprae organisms in mice treated by BDQ administered orally

Table 1 .

 1 Multiplication of M. leprae organisms in nude mice to determine minimal effective dose of BDQ administered orally active against M. leprae and the benefit of the combination of CFZ and BDQ.

	Treatment a	Positive footpad/	Range of AFB/positive footpads	p value b
		total footpad	(mean log 10 AFB/ footpad
			standard deviation on positive
			footpads)
	Untreated control	10/10	7.79-8.67 (8.13±0.26)	/
	RIF 10 mg/kg	6/10	4.54-5.72 (NA c )	0.08
	BDQ 0.10 mg/kg	10/10	6.62-8.06 (7.39±0.50)	1
	BDQ 0.33 mg/kg	10/10	6.50-7.57 (7.02±0.40)	1
	BDQ 1 mg/kg	10/10	7.08-7.94 (7.56±0.32)	1
	BDQ 3.3 mg/kg	0/10	NA c	0.00001
	BDQ 25 mg/kg	0/10	NA c	0.00001
	CFZ 20 mg/kg	3/10	4.54 (NA c )	0.003
	BDQ 0.10 mg/kg + CFZ 20 mg/kg	2/10	4.54-5.15 (NA c )	0.0007
	BDQ 0.33 mg/kg + CFZ 20 mg/kg	1/10	5.54 (NA c )	0.00009
	a treatment was given by oral gavage 5 days a week for 24 weeks beginning one month after inoculation,

Table 2 .

 2 Bactericidal activity against M. leprae THAI53 of bedaquiline administered orally and bedaquiline long-acting measured in Swiss mice by the proportional bactericidal method.

	Treatment a	No. of footpads showing multiplication b of M. leprae/ No. of footpads harvested, by inoculum 5x10 3 5x10 2 5x10 1 5x10 0	% viable M. leprae c	% viable M. leprae killed by treatment d	p value e	p value h
	Untreated control	10/10	8/10	8/10	2/10	2.753	/
	RIF 10 mg/kg SD f RIF 10 mg/kg 3D g	8/10 0/10	6/10 0/10	4/10 0/10	--	0.275 0.004	90.01 99.85	0.0046 <0.001 <0.001
	BDQ 2 mg/kg SD BDQ-LA IM 2 mg/kg SD	9/10 9/10	7/10 8/10	6/10 5/10	--	0.692 0.692	74.86 74.86	0.075 0.069	1
	BDQ 2 mg/kg 3D BDQ-LA IM 2 mg/kg 3D	9/10 8/10	5/10 6/10	5/10 3/10	--	0.347 0.219	87.39 92.04	0.009 0.0016	0.584
	BDQ 20 mg/kg SD BDQ-LA IM 20 mg/kg SD	3/10 7/10	0/10 3/10	0/10 1/10	--	0.009 0.055	99.67 98.00	<0.001 <0.001	0.005
	BDQ 20 mg/kg 3D BDQ-LA IM 20 mg/kg 3D	0/10 3/10	0/10 2/10	0/10 0/10	--	0.004 0.014	99.85 99.49	0.001 <0.001	0.014

b comparisons of the proportion of mice with positive footpads of each treated group versus untreated control (a p-value <0.005 was considered to be statistically significant when applying Bonferroni's correction).. c not applicable

Author summary

The current multidrug therapy is effective against leprosy but remains long and difficult to observe for patients supporting the need of monthly -based treatment. Bedaquiline (BDQ), a diarylquinoline with a long half-life, is a candidate drug to shorten leprosy treatment by targeting the electron transport chain and inhibiting the ATP synthesis. In this work, we demonstrated that (i) the minimal effective dose of BDQ administered orally against M. leprae is 3.3 mg/kg, (ii) BDQ did not add significantly to the efficacy of CFZ at the doses tested, and (iii) BDQ long acting formulation was similar or less active than BDQ administered orally at equivalent dosing and frequency but should be tested at higher dosing in order to reach equivalent exposure in further experiments.

Discussion

The current length of the leprosy treatment remains a challenge [START_REF] Lockwood | Chapter 2.6 Treatment of leprosy[END_REF]. During the last decades, few new antituberculous drugs were synthesized and even rarer are those active against M. leprae. In 2006, Ji et al. demonstrated the antileprosy activity of BDQ administered orally, a new diarylquinoline which was also active against other mycobacteria such as M. tuberculosis [START_REF] Ji | Bactericidal activities of R207910 and other newer antimicrobial agents against Mycobacterium leprae in mice[END_REF][START_REF] Andries | A diarylquinoline drug active on the ATP synthase of Mycobacterium tuberculosis[END_REF]. They showed that a single dose of 25 mg/kg BDQ was as effective as rifapentine, moxifloxacin, as well as RIF, which is currently the most powerful antileprosy drug. Nevertheless, the minimal effective dose of the BDQ administered orally is unknown despite knowing that it may enable developing a safe dose. In our present work, we determined that a dose as low as 3.3 mg/kg of BDQ administered orally displayed a bactericidal activity indistinguishable from a 25 mg/kg dose (Table 1, Figure 1). The MED of 3.3 mg/kg is higher than the ≤1 mg/kg that was found to be the lowest active dose in an immunocompetent mouse model of leprosy [START_REF] Gelber | The Diarylquinoline R207910 Is Bactericidal against Mycobacterium leprae in Mice at Low Dose and Administered Intermittently[END_REF] which can be explained by the lack of T-cell immunity, known to be important in mycobacterial diseases, in the nude mouse model used to determine the MED of BDQ administered orally in the present study. The MED of BDQ administered orally was not achieved in the immunocompetent mouse model because all the doses tested including the lowest one (1 mg/kg) were able to reach the limit of detection of M. leprae present in mouse footpads as enumerated by the AFB microscopy [START_REF] Gelber | The Diarylquinoline R207910 Is Bactericidal against Mycobacterium leprae in Mice at Low Dose and Administered Intermittently[END_REF] suggesting a MED of ≤ 1mg/kg. The MED of 3.3 mg/kg obtained in the immunocompromised mouse model of leprosy is lower than that obtained with M. tuberculosis in an immunocompetent mouse model (6.5 mg/kg) [START_REF] Gelber | The Diarylquinoline R207910 Is Bactericidal against Mycobacterium leprae in Mice at Low Dose and Administered Intermittently[END_REF] but the MED of BDQ administered orally against M. tuberculosis should be compared to the MED of BDQ administered orally against M. leprae determined in an immunocompetent mouse model (MED ≤ 1 mg/kg) suggesting at least a 6-


Sébastien Guillet
email: sebastien.guillet@dendrolab.ch

Christophe Corona

Markus Stoffel
email: markus.stoffel@dendrolab.ch

Myriam Khodri

Franck Lavigne

Pablo Ortega

Nicolas Eckert

Valérie Pascal Dkengne Sielenou

Sébastien, Corona, Christophe, Stoffel Markus Guillet

Franck, Ortega Pablo Lavigne

Eckert Nicolas Sielenou

Pascal Dkengne

Valérie Daux

Olga V Churakova

Nicole Davi

Jean-Louis Edouard

Yong Zhang

Brian H Luckman

Vladimir S Myglan

Joël Guiot

Martin Beniston

Valérie Masson-Delmotte

Clive Oppenheimer

Climate response to the Samalas volcanic eruption in 1257 revealed by proxy records

historical evidence exists for persistent dust veils and widespread NH climatic cooling [START_REF] Büntgen | Cooling and societal change during the Late Antique Little Ice Age from 536 to around 660 AD[END_REF][START_REF] Briffa | Influence of volcanic eruptions on Northern Hemisphere summer temperature over the past 600 years[END_REF][START_REF] Stoffel | Estimates of volcanic-induced cooling in the Northern Hemisphere over the past 1,500 years[END_REF] . This absence of evidence of strong cooling in 1258 even called into question the ability of tree-ring width proxies to detect the short-term cooling [START_REF] Mann | Underestimation of volcanic cooling in tree-ring-based reconstructions of hemispheric temperatures[END_REF][START_REF] Mann | Discrepancies between the modeled and proxyreconstructed response to volcanic forcing over the past millennium: Implications and possible mechanisms[END_REF][START_REF] Anchukaitis | Tree rings and volcanic cooling[END_REF][START_REF] D'arrigo | Volcanic cooling signal in tree ring temperature records for the past millennium[END_REF] associated with the largest volcanic eruptions and the skill of climate models to predict the magnitude and persistence of volcanic cooling if they do not incorporate aerosol microphysical modules [START_REF] Stoffel | Estimates of volcanic-induced cooling in the Northern Hemisphere over the past 1,500 years[END_REF][START_REF] Timmreck | Limited temperature response to the very large AD 1258 volcanic eruption[END_REF] .

Historical evidence for a strong cooling in 1258

Here, we shed light on the climate forcing associated with the 1257 Samalas eruption through the analysis of an extensive compilation of mediaeval texts, drawn from the Monumenta Germaniae Historica, the Rerum Britannicarum Medii Aevi Scriptores and the Recueil des historiens des Gaules et de la France (Text S1). A total of 35 narrative sources (Table S1) attest to significant and widespread climate anomalies across Western Europe in the spring, summer and autumn of 1258 (Fig. 1). Notably, such abundant and detailed documentation in mediaeval archives is reserved for extreme meteorological events: the only comparable episode, reported in more than 30 sources (Text S1), relates to an extremely cold winter in Europe in 1233/34 [START_REF] Pfister | Winter air temperature variations in western Europe during the Early and High Middle Ages (AD 750-1300)[END_REF] .

Without exception, all chronicles surveyed report cold, incessant rainfall and unusually high cloudiness in 1258, which prevented crops and fruits from reaching maturity. The Norman Notes of Coutances provides a good example: "There was no summer during summer. The weather was very rainy and cold at harvest time, neither the crop harvest nor the grape harvest were good. Grapes could not reach maturity; they were green, altered and in poor health". This prefigures abundant descriptions of impacts on crops of the "Year Without a Summer" in 1816, a year after the Tambora eruption [START_REF] Oppenheimer | Climatic, environmental and human consequences of the largest known historic eruption: Tambora volcano (Indonesia) 1815[END_REF][START_REF] Stommel | Volcano weather: the story of 1816, the year without a summer[END_REF][START_REF]The Year without a summer?: world climate in 1816[END_REF] S1.

Holy Roman Empire

Poland

Kingdom of France

Sweden

Norway

Serbia

Kingdom of Leon and Castille

Kingdom of England

Kingdom of Hungaria

Kingdom of Aragon

Kingdom of Portugal Kingdom of Sicily

Papal States

Kingdom of Denmark

Emirate of Grenade

State of the Teutonic Order

Province of Guyenne

Iceland

Ireland

Kingdom of Scotland

Prague (Czech Republic), but quantitative information on grape harvest dates (GHD) could only be retrieved for France. In Alsace-Lorraine, Richer of Senones reports that grapes were still green and as "hard as stone" on 8 October 1258 and that the harvest was postponed until the end of October (i.e., circa day 305 of the year; doy), as winemakers were hoping for improved weather conditions in autumn (Methods and Text S2). In Paris, the Annales Clerici Parisiensis point to a GHD start on doy 312 (Fig. 2a, Fig 3b), whereas in Burgundy, the start of GHD was estimated to doy 302 (Methods and Texts S2, S3). If compared with existing GHD series from Alsace-Lorraine (1700-2005), Ile-de-France (1478-1977) and Burgundy (1354-2006), we find that these are by far the latest dates ever observed, and that grape harvest started 7 (doy 298), 16 (doy 312), and 6 (doy 296) days later than in the 1816 "year without a summer". Using Extreme Value Theory, we find that the return period of the 1258 GHD exceeds 1000 yr for the three regions (as compared to circa 100 yr for 1816) and that it reaches infinite values if known volcanically-influenced years are excluded from the analyses. This suggests that their occurrence is exceedingly unlikely to be related to interannual climate variability (Methods and Text S3). We caution, however, that the GHD approach does not account for possible changes in grape varieties or in agricultural and/or anthropogenic practices over time. As such, the return periods here are merely indicative of the unusual character of the grape growing season -which is controlled by maximum April-September temperatures [START_REF] Daux | An open-database of Grape Harvest dates for climate research: data description and quality assessment[END_REF] -over Western Europe in 1258 (Fig. 2b). The likely volcanic origin of the extreme cooling is reinforced by two accounts attesting the presence of a persistent dust veil (Text S4) over Europe in 1258. In Germany, the Annals of Speyer reveal that this year was commonly referred to as a munkeliar (Fig. 3a), meaning dark year or year of fog in contemporary German. The presence of a dense dust veil is corroborated by a detailed testament of a very dark (Danjon L=0; Methods, Text S4), total lunar eclipse [START_REF] Keen | Volcanic Aerosols and Lunar Eclipses[END_REF]

Magnitude of NH summer cooling inferred from proxy records

Complementary information is provided by paleoclimate archives. We have used a tree-ring network covering the period 1000-2000 CE to quantify and map climatic anomalies induced by the Samalas eruption and, in particular, for those regions of the NH which remain unrepresented by mediaeval archives (Fig. S1, Table S2). All of the 25 chronologies span the 13 th century and are significantly correlated with summer (JJA) temperatures (Fig. S2). A JJA NH temperature reconstruction (40-90°N over land) was generated using a nested approach containing 12 nests that account for the decrease in the number of regional chronologies back in time. Each nest passed all verification tests (R 2 and RE for the calibration period, as well as r 2 and CE for the verification period, see Methods and Table S3) for the period 1805-1972. The most replicated nest (1230-1972, n=25) accounts for 41% of the variance in the instrumental data with RE=0.37 and CE=0.34, demonstrating the robustness of the reconstruction for the 13 th century. To quantify cooling with respect to contemporary climatology, the reconstruction was filtered with a 30-yr running mean (see Methods). We find that the cumulative distribution functions for all major volcanic eruptions and non-volcanic years differ in a statistically significant manner (Kolmogorov-Smirnov test, p<0.001) (Fig. 4a). Extreme average surface cooling occurred in 1259, reaching -1.2°C, making it the third coldest summer since at least 1200 CE, whereas 1258 (-0.7°C), 1601 (-1.2°c), 1453 (-1°C), and 1816 (-1.2°C) rank as the 16 th , 2 nd , 7 th , and 1 st coldest years in the series, respectively. If JJA temperature anomalies are summed for the two years after an eruption The same year, wine, wheat and other fruits were greatly altered and this year was also commonly referred to as munkeliar". The use of the Middle High German expression munkeliar, rather than its Latin equivalent (annus obscuritatis or annus caliginis), suggests that the exceptional persistence and intensity of insolation dimming was not only omnipresent but unusual enough for commoners to give it a proper name (Source: Speyrer Kopialbuch. Generallandesarchiv Karlsruhe GLA 67 Nr. 448 fol. 39v). b, Contemporary illustration of wine harvesting as illustrated in the Martyrology of the Saint-Germain-des-Prés Abbey (Source: National Library of France, Paris, Ms lat. 12834, fol. 69v).

(Fig. 4b), the Samalas eruption is seen to be associated with the strongest cooling anomaly of the last millennium.

Spatial heterogeneity of NH summer cooling

The 25 chronologies were grouped into 11 clusters to estimate regional variability of summer cooling induced by the Samalas eruption, and complemented with three annually-resolved, summer season stable oxygen isotope δ 18 O series from Greenland ice cores (GRIP, Crete, DYE3). The spatial extent of each cluster was defined via spatial fields of correlation between each proxy record and mean JJA temperatures of the Berkeley Earth Surface Temperature (BEST) dataset [START_REF] Rohde | A New Estimate of the Average Earth Surface Land Temperature Spanning 1753 to 2011[END_REF] (Fig. S2). Based on this approach, Greenland ice cores provide information on temperature in Eastern Greenland and Iceland.

The JJA gridded temperature reconstruction (3486 1°×1° grid points with a CE>0.1, Fig. S3) is given for 1258, 1259, 1453, 1601, and 1816 (Fig. 4c).

For all these years, volcanic forcing induced pronounced cooling over Siberia and, to a lesser extent, in Western Europe. Less consistent patterns are observed in Alaska and on the Western Pacific Coast. In 1258, Western Europe, the Canadian Rockies, Central Asia and Siberia were the most severely affected with a reconstructed cooling between -1.4 and -2°C (wrt 30-yr climatology). The extremely cold conditions prevailing in Siberia and Central Asia are further corroborated by the presence of frost rings, which have previously been recognized as a response to severe cooling associated with volcanic forcing of climate [START_REF] Lamarche | Frost rings in trees as records of major volcanic eruptions[END_REF] (Table S4).

In contrast, the Greenland ice core cluster (-0.1°C, rank 408) and tree-ring records from Quebec (+0.5°C, rank 676), Alaska (+0.3°C, rank 617) and, to a lesser extent, Scandinavia (-0.5°C, rank 264) do not reveal pronounced cooling, suggesting that, in these regions, internal modes of climate variability outweighed the direct radiative effects of volcanic aerosol. The warm anomalies over Alaska in 1258-9 (Fig. 4c) could arise from a positive phase of the Pacific-North-Atlantic pattern, linked to ENSO (Fig. S4). Climate simulations and proxy records indeed show increased probability for an El Niño event to occur in the first or second year after a large volcanic eruption [START_REF] Emile-Geay | Volcanoes and ENSO over the Past Millennium[END_REF][START_REF] Adams | Proxy evidence for an El Niño-like response to volcanic forcing[END_REF][START_REF] Li | Interdecadal modulation of El Niño amplitude during the past millennium[END_REF] . This assumption is further confirmed by El Niño conditions inferred for 1258/59 from tree-ring and sediment proxies [START_REF] Emile-Geay | Volcanoes and ENSO over the Past Millennium[END_REF] .

We also use tree-ring records and historical archives to investigate the persistence of volcanicallyinduced cooling. Interestingly, Fig. 4c shows that several regions experienced maximum cooling only in 1259, i.e., two years after the eruption. In Quebec, for instance, the strong negative anomaly in 1259 likely results from extremely cold conditions at the end of the previous growing season, as revealed by the presence of light rings in the latewood of 1258 (Table S4). This lag effect, known as autocorrelated biological memory [START_REF] D'arrigo | Volcanic cooling signal in tree ring temperature records for the past millennium[END_REF] , is often observed in tree-ring width (TRW) series and reflects the influence of climatic and physiological conditions prevailing in the year prior to tree-ring formation (e.g., via needle generation, or the creation of carbohydrate reserves). In Siberia, however, biological memory cannot fully explain the temperature drop in 1259, because frost rings in the Yamal, Polar Ural, Altai, and Mongolian tree-ring chronologies (Table S4) as well as historical records documenting abundant snowfall in the Altai in July clearly confirm the extreme cooling experienced in these regions in summer 1259 (ref. 33).

Based on the available data, the only regions to experience mild conditions in the summer of 1259 were the Pacific coast of North America (as identified in tree-ring records) and Western Europe, for which five historical sources were found (Text S6). In the Annales Wormatienses (Germany), it is stated that on the day of the Apostles Peter and Paul [29 June], the summer was hot and dry. From March to August, little or even no rain fell, however wine and all other fruits of the earth were abundant. Treering proxies and chronicles indicate a general reduction in all climatic anomalies by 1260 and 1261 (Text S7), and thus provide strong evidence that the impacts of the eruption ceased after 3-4 yr (Fig. S5). This points to a more ephemeral effect on climate than suggested by computational models that exclude aerosol microphysics, which predicted a persistence of temperature anomalies over extratropical land of the order of -1°C until 1264 (ref. 34).

The estimated hemispheric cooling following the Salamas eruption is comparable to that observed for 1453, 1601, and 1816, although Samalas released substantially more sulphur to the atmosphere than the other volcanic events. This underlines the hypothesis that volcanic cooling does not increase linearly with the total sulphur yield of an eruption [START_REF] Timmreck | Limited temperature response to the very large AD 1258 volcanic eruption[END_REF][START_REF] Pinto | Self-limiting physical and chemical effects in volcanic eruption clouds[END_REF] . In addition, we confirm that the magnitude and persistence of cooling resulting from the 1257 event were much smaller and shorter, respectively, than predicted in PMIP3 simulations. The disparity is sufficiently great that the Samalas eruption can no longer be held responsible for a dramatic and persistent cooling of continental regions of the NH [START_REF] Mann | Underestimation of volcanic cooling in tree-ring-based reconstructions of hemispheric temperatures[END_REF] .

Evidently, the largest eruption of the Common Era only induced "years without summer" in certain regions, but not across the entire NH. This finding is consistent with observations in 1816, and points to the modulation of radiative forcing effects of volcanic aerosols by internal climate variability (e.g., ENSO, NAO) [START_REF]The Year without a summer?: world climate in 1816[END_REF] . Our reconstructed patterns of spatially heterogeneous temperature anomalies provide crucial constraints for the evaluation of climate model performance for the 1257 Samalas event and other high sulphur yield eruptions in the tropics.

Complex societal response to the Samalas eruption

It has been posited that a cluster of volcanic eruptions in 536, 540 and 547 may have played a crucial role in the rise and fall of several polities, pandemics and human migration [START_REF] Büntgen | Cooling and societal change during the Late Antique Little Ice Age from 536 to around 660 AD[END_REF] . Several sources point to severe food shortages and subsistence crises in parts of Europe (Kingdom of France, Kingdom of England, Holy Roman Empire, Iberian Peninsula) in 1258 and 1259, although they were less widespread and shorter than the famines of 1195-1197, 1233-1235 or 1315-1319 (ref. 36,37), which cannot be associated with volcanic events. The most severe socio-economic consequences reported at the time of the Samalas eruption are from England, where the 1258 famine caused by two consecutive years of bad harvests (1256-57), high prices and speculation may have killed about 15,000 people in London alone [START_REF]A bioarchaeological study of medieval burials on the site of St Mary Spital: excavations at Spitalfields Market[END_REF] . Likewise, the Shôga famine (Japan, 1257-1260) was possibly amplified by adverse climatic conditions in 1258 and 1259 (ref. 26). In both cases, we conclude that volcanic forcing of the climate following the Samalas eruption played a role in aggravating these crises, but should not be considered as the trigger for the famines or related social unrest. This emphasizes the complexity of identifying causal connections between individual cases of volcanically induced climate anomalies and major societal upheavals [START_REF] Oppenheimer | Eruption politics[END_REF] .

Methods

Analysis of historical sources. We exhaustively reexamined 180 annals and chronicles (original and edited manuscripts, Table S1, Figs. 12), written in the 13 th century to document the climatic impacts of the Samalas eruption. Sources were mostly in Latin and to a lesser extent in the vernacular (Old French, High Middle German, and Old Castillan). From this corpus, 59 narrative sources were from Western Europe and report weather and phenological (crop and grape harvest dates) observations for the years 1257, 1258, 1259, 1260 and 1261 (Table S1). Particular attention was also paid to historical reports of optical phenomena attesting to the presence of volcanic aerosols in the stratosphere [START_REF] Keen | Volcanic Aerosols and Lunar Eclipses[END_REF] (e.g., dimming of the Sun, dark total lunar eclipse, Text S4). The reliability of eclipse observations was further assessed through the catalogue of five millennia of lunar eclipses (1999 BCE -3000 CE) [START_REF] Espenak | Five Millennium Catalog of Lunar Eclipses[END_REF] . The darkness of the moon was rated according to the Danjon scale [START_REF] Danjon | Relation Entre l'Eclairement de la Lune Eclipsée et l'Activité Solaire[END_REF] ranging from L=0 (very dark) to L=4 (very bright copper-red or orange). The focus was clearly on first-hand information derived from contemporary sources, and dates in chronicles were screened systematically for errors by cross-checking information with other manuscripts.

Analysis of grape harvest dates. Mediaeval grape harvest dates (GHD) obtained from Saint-Denis and Senones

Abbeys were added to GHD series available for the French regions of Alsace and Ile de France [START_REF] Daux | An open-database of Grape Harvest dates for climate research: data description and quality assessment[END_REF] . Mediaeval dates were first converted to Gregorian calendar by adding 7 days to Julian calendar and then expressed as day of the year (doy). Correlation analyses were performed against April-September temperatures from Paris-Montsouris (1676-1977) (ref. 42), Dijon (1951-2006) and Strasbourg (1950-2005) (https://climexp.knmi.nl/) to assess the suitability of GHD as a temperature proxy (Fig. 2). Missing values were calculated on the basis of inter-series similarity, allowing estimation of the 1258 GHD in Burgundy. For all three series, extreme GHD were selected as the 10% latest dates. The three estimated datasets as well as a composite record were fitted using the discrete version of the Generalized Pareto Distribution (GPD), leading to four estimates of the 1258 GHD return period. Also, the analysis was duplicated, excluding the sample of years corresponding to known large volcanic eruptions from the GPD fitting, providing estimates of the 1258 GHD return period under unperturbed climate conditions (Text S3).

Tree-ring chronology selection and assessment. We compiled an extensive database of tree-ring width (TRW and maximum latewood density (MXD) chronologies from the International Tree-Ring Database (ITRDB) and published papers from NH tree-line sites located between latitudes 39 and 73°N (Table S5). For each individual chronology, we carefully tested its sensitivity to June-August (JJA) temperatures (Fig. S2) and excluded series with an ambiguous climate response. The remaining 25 chronologies (13 TRW and 12 MXD) all include the 13 th century; they were grouped into 11 regional clusters based on correlation matrices (r>0.3, p-value<0.05). All clusters document the Samalas eruption with an expressed population signal or EPS never <0.85.

Tree-ring based NH temperature reconstruction. Each cluster was standardized using standard (i.e., negative exponential), RCS [START_REF] Esper | Low-Frequency Signals in Long Tree-Ring Chronologies for Reconstructing Past Temperature Variability[END_REF][START_REF] Helama | Regional curve standardization: State of the art[END_REF] and Signal-Free RCS [START_REF] Melvin | A 'signal-free' approach to dendroclimatic standardisation[END_REF] standardizations that preserve low-frequency (that is, multi-decadal) temperature changes. We then transferred this record into JJA temperature anomalies (wrt 1961-1990) through a bootstrap linear model using principal component analysis calibrated against JJA land surface temperature anomalies (1805-1972) from the BEST dataset [START_REF] Rohde | A New Estimate of the Average Earth Surface Land Temperature Spanning 1753 to 2011[END_REF] . The calibration and validation process was repeated 1,000 times using a bootstrap method to assess the robustness of the transfer function. To test the quality of the fit of each individual run, we used coefficient of determination statistics (R 2 for the calibration and r 2 for the verification periods), whereas the RE (reduction of error) and CE (coefficient of efficiency) statistics served to test the predictive capacity of the transfer function. Calibration and validation statistics are illustrated for each nest with their 2.5 and 97.5 percentiles (Table S3) and the reconstruction is given with its 95%-confidence intervals (Text S5). To account for the decreased number of chronologies back in time, we used a nested approach (12 nests) to obtain the longest possible reconstruction. The final reconstruction was developed by splicing all the nested time series after adjustment of mean and variance of each nested reconstruction segment to the best replicated nest (1230-1972). This approach is known to stabilize the variance of the final reconstruction but to increase the error variance. Calibration and validation statistics (R 2 =0.39-0.45, r 2 =0.35-0.38; reduction of error RE=0.36-0.39, coefficient of efficiency CE=0.34-0.37; see Table S3) were used to evaluate the reliability of the reconstruction. In Fig. 4b, reconstructed temperature anomalies are presented as deviations from the climate mean climatology smoothed with a 15-yr time window on either side of the volcanic eruptions dated to 1257, 1452, 1600, and 1815.

For example, in the case of Samalas, a background was calculated by averaging the window 1243-1257 and 1259-1273. The anomaly is then created by subtracting this background from the 1258 reconstructed temperature. As for GHDs, the 10% lowest values were fitted using the continuous GPD to evaluate the return period of the cooling induced by the Samalas and by other big eruptions. Analysis was done for the 1258 and 1259 annual values, as well as for anomalies cumulated over two years (Fig. 4).

Regional reconstructions.

To account for the spatial variability of volcanic cooling, the 25 tree-ring chronologies were grouped into 11 clusters. In addition, to extend the analysis to other unrepresented regions we included a cluster of three stable isotope series from Greenland ice cores (GRIP, Crete, DYE3), all statistically significantly related to summer temperature conditions over Eastern Greenland and Iceland [START_REF] Vinther | Climatic signals in multiple highly resolved stable isotope records from Greenland[END_REF] . The spatial extent of each cluster was defined via the spatial field of correlations between each proxy record and mean JJA temperatures of the BEST dataset (Fig. S2). We then applied linear regression analysis (for clusters with only one chronology) or a bootstrap linear model using principal component analysis (for clusters with multiple chronologies) in order to calibrate cluster series to JJA gridded temperature anomalies (wrt 1961-1990) from the BEST dataset. Full period calibration was performed separately for each of the 12 clusters over the period 1901-1990 while the robustness of this association was tested through cross-validation procedures for the period 1901-1945 and 1946-1990 (Fig. S3).

Only grid points with a CE>0.1 were used for the spatial representation of temperature anomalies.

Figure 1 .

 1 Figure 1. Spatial extent of weather and optical anomalies observed in Europe in 1258. All sources are listed in TableS1.

Figure 2 .

 2 Figure 2. Grape harvest dates in France (1258-2006). a, Continuous records of days of year (doy) on which GHD occurred for Ile-de-France (purple), Alsace (orange) and Burgundy (green) between 1350 and 2006. Data are complemented by newly discovered sources for the years 1258, 1279, and 1294. The green star indicates the estimated 1258 GHD in Burgundy (Text S3). The most delayed grape harvest of the last 800 years occurred in 1258. See Table S1 for a list of all sources. b, GHD for Ile-de-France, Alsace, and Burgundy are significantly correlated with April-September mean air temperatures of the Paris-Montsouris, Strasbourg and Dijon meteorological stations.

  on 12 November 1258 in the Annales Ianuenses (Genoa, Italy). It complements the previously known dark lunar eclipse described in the Chronicle of John de Taxster (England) for 23 May 1258 (ref. 5). Beyond Europe, only a very limited number of contemporary narrative sources have hitherto been surveyed for evidence of climate anomalies. The Mirror of the East (Azuma Kagami, Kamakura, Japan) is an exception 26 and reports a wet and cold summer accompanied by heavy rain and strong winds, which reportedly destroyed paddies and gardens in several provinces.

Figure 3 .

 3 Figure 3. Original contemporary manuscript from the annals of Speyer describing the dust veil and climate anomalies observed in 1258 CE. a, The text in Latin says: "The same year, wine, wheat and other fruits were greatly altered and this year was also commonly referred to as munkeliar". The use of the Middle High German expression munkeliar, rather than its Latin equivalent (annus obscuritatis or annus caliginis), suggests that the exceptional persistence and intensity of insolation dimming was not only omnipresent but unusual enough for commoners to give it a proper name (Source: Speyrer Kopialbuch. Generallandesarchiv Karlsruhe GLA 67 Nr. 448 fol. 39v). b, Contemporary illustration of wine harvesting as illustrated in the Martyrology of the Saint-Germain-des-Prés Abbey (Source: National Library of France, Paris, Ms lat. 12834, fol. 69v).

Figure 4 .

 4 Figure 4. Tree-ring reconstructions of NH extra-tropical land (40-90° N) summer temperature anomalies since 1000 CE. a, Summer (JJA) temperature anomalies following the 1257 Samalas eruption in 1258 (red) and 1259 (blue) as compared to cumulative distribution functions (cdf) for all major volcanic eruptions (i.e., 1109, 1453, 1601, 1641, 1695, 1783, 1809, 1816, 1835, 1884, and 1912; black) and for all non-volcanic years (red) since 1000 CE. b, same as in a, but for groups of two consecutive years following major eruptions. c, Spatial extent of the JJA temperature anomalies induced by the Samalas (cooling shown for 1258 and 1259), unknown (1453), Huaynaputina (1601), and Tambora (1816) eruptions. For details see Methods.

  wrt 31yr mov. av.)

Acknowledgments

S.G., C.C., M.S. and O.C.S. acknowledge support from the Era.Net RUSplus project ELVECS (SNF project number: IZRPZ0_164735). This study benefited from data gathered within the ANR CEPS GREENLAND project.

V.S.M. received support from the Russian Science Foundation (project № 15-14-30011). Rashit Hantemirov kindly provided a millennium-long chronology. The authors are grateful to William S. Atwell and William Wayne-Farris for discussions on historical sources from Japan as well as to Marie Luisa Avila for her help with Muslim sources from Mediaeval Spain. S.G. and C.C. are very grateful to Séverine Finet, Laurence Fazan and Pascal Guérin for their help with R-scripts, translations and fruitful discussions, respectively. The authors are very grateful to the referees and senior editor Dr. Amy Whitchurch for their insightful comments, which greatly improved the manuscript.

Data availability. The data used to perform our analysis as well as our results can be accessed at: https://www.ncdc.noaa.gov/data-access/paleoclimatology-data/datasets/climate-reconstruction

Author contributions

Additional information

Supplementary information is available in the online version of the paper Correspondence and requests for materials should be addressed to S.G.

Competing financial interests

The authors declare no competing financial interests.


A Urdaneta-Calzadilla

N Galopin

I Niyonzima

Olivier Chadebec

G Meunier

B Bannwarth

FEM-BEM modeling of nonlinear magnetoelectric effects in heterogeneous composite structures

Keywords: Boundary element methods, Finite element analysis, Magnetoelectric effects, Rotating coilless ME device

published or not. The documents may come

FEM-BEM modeling of nonlinear magnetoelectric effects in heterogeneous composite structures

I. INTRODUCTION

E LECTROMAGNETIC interaction can be achieved by magnetoelectric (ME) composites. These materials are based on mechanical exchanges through the combination of piezoelectric and magnetostrictive materials which have, respectively, strong electro-mechanical and magneto-mechanical couplings [START_REF] Eerenstein | Multiferroic and magnetoelectric materials[END_REF]. Numerical modeling of these phenomena is usually performed through the Finite Element Method (FEM) [START_REF] Nguyen | Finite element modeling of magnetic field sensors based on nonlinear magnetoelectric effect[END_REF], [START_REF] Talleb | Finite element modeling of magnetoelectric laminate composites in considering nonlinear and load effects for energy harvesting[END_REF]. Nevertheless, to correctly account for the decay of magnetic fields at infinity, the FEM needs a big enough free space domain to be explicitly considered and meshed. As a result, a large number of degrees of freedom being located outside the active material, notably with regard to the quantity of active materials and their distance from the sources. This can be avoided by setting up a coupling between the Boundary Element Method (BEM) and the FEM for the magnetic part, and the FEM for the electric and mechanical parts. The coupling between FEM and BEM has proved to be a powerful tool for reducing computational time and allocated memory for solving multi-physics problems, and enabling the use of a single mesh for all the problems concerned (see e.g., [START_REF] Urdaneta-Calzadilla | A FEM-BEM coupling strategy for the modeling of magnetoelectric effects in composite structures[END_REF] by the authors). It should be noted, however, that this previous work was limited to a linearized approach to magneto-mechanical effects, which severely limited the validity of the model.

Contrary to piezoelectricity which in most situations can be accurately described by linear relations, magneto-mechanical behavior is strongly nonlinear. Many approaches have been taken to consider this non-linearity. In this paper, we implemented the description of magneto-mechanical behavior in terms of a scalar invariant's formulation of the Helmholtz free energy [START_REF] Ahmed | 3-D Magneto-Mechanical Finite Element Analysis of Galfenol-Based Energy Harvester Using an Equivalent Stress Model[END_REF]. This approach has already been used in a purely FEM context. The main contribution of this article is to extend it to FEM-BEM coupling. We will see that this extension leads to numerical difficulties, which have been resolved by using an adaptive iterative approach of nonlinear single-physics problems.

In the rest of the paper, we present the material laws and the formulations used in the model. Then the model is used to simulate a rotating coilless ME device under the influence of a permanent magnet. This device can operate as an energy harvester or as an actuator. In particular, in actuator mode, the calculation of magnetic forces is discussed using two approaches: the magnetic charge method and the Laplace force.

II. CONSTITUTIVE LAWS

The considered behavioral law for the piezoelectric phase is given by the two linear coupled relationships (1),

T (S, E) = C E : S -t e • E, D(S, E) = e : S + ε S • E, (1)
with D the electric displacement field, E the electric field, T the Cauchy stress tensor, S the linear strain tensor, C E the elasticity tensor at constant electric field, e the piezoelectric tensor and ε S the permittivity tensor at constant strain." • ", " : " and t • denote respectively the dot product, the double dot product and the transpose operator. Inside the magnetostrictive phase, the Helmholtz free energy ψ can be expressed in the case of isotropic materials as a polynomial expression on its invariants I i [START_REF] Ahmed | 3-D Magneto-Mechanical Finite Element Analysis of Galfenol-Based Energy Harvester Using an Equivalent Stress Model[END_REF]:

ψ = 1 2 λI 2 1 + µI 2 + nα i=1 α i I i 4 + n β i=1 β i I i 5 + nγ i=1 γ i I i 6 , (2)
where:

I 1 = tr(S), I 2 = tr(S 2 ), I 4 = B • B I 5 = B • S • B, I 6 = B • S2 • B, (3)
with B the magnetic induction, H the magnetic field, S the deviatoric part of the strain tensor, λ and µ the Lamé coefficients of the material, α i , β i and γ i a finite set of scalar coefficients. The magnetic and mechanical coupled nonlinear behavioral laws can then be obtained by the analytical differentiation of the energy ψ:

T (B, S) = ∂ψ ∂S , H(B, S) = ∂ψ ∂B . ( 4
)
The behavior of ME composite structures is then obtained through the combination of (1) and (4).

III. FORMULATIONS AND RESOLUTION METHOD

The electrical and mechanical problems are limited to the active domain Ω m [START_REF] Urdaneta-Calzadilla | A FEM-BEM coupling strategy for the modeling of magnetoelectric effects in composite structures[END_REF]. Indeed, because of the high permittivity of piezoelectric materials and since the electrodes are in direct contact with the piezoelectric materials, the electric displacement field is canalyzed inside the materials. Therefore, restricting the domains to Ω m is a reasonable approximation. Considering static behaviour, from Maxwell-Gauss relation, [START_REF] Eerenstein | Multiferroic and magnetoelectric materials[END_REF] and appropriate boundary conditions, the electric weak form in terms of electric scalar potential (E = -∇φ) reads, find φ such that:

Ωm ∇δφ • D(S, E) dΩ = 0 ∀δφ. ( 5
)
with δφ an appropriate virtual electric scalar potential field.

From ( 1), ( 4) and the balance of linear momentum, the mechanical weak form in terms of mechanical displacement u reads, find u such that:

Ωm ∇ S δu : T (S, E, B) dΩ = 0 ∀δu, (6)
where ∇ S u = 1 2 (∇u + (∇u) t ) and δu an appropriate virtual mechanical displacement vector field. Both ( 5) and ( 6) are discretized using first order FEM and nodal shape functions. Similarly to [START_REF] Urdaneta-Calzadilla | Modeling of Magnetoelectric Effects in Composite Structures by FEM-BEM Coupling[END_REF], the weak form of the magnetic problem is obtained by testing Maxwell-Ampere's equation, keeping the boundary term. It reads, find (a, ϕ red ) such that:

Ωm ∇ × δa • H(B, S) dΩ + ∂Ωm ∇ × δa • n ϕ red d∂Ω = ∂Ωm (δa × n) • H 0 d∂Ω ∀δa, (7)
with a the magnetic vector potential (B = ∇ × a), δa an appropriate virtual magnetic vector potential field , ϕ red the magnetic reduced scalar potential, H 0 the magnetic field created by the currents external to the domain calculated by the Biot-Savart law and n the outwards normal vector to ∂Ω m . The previous equation is discretized using edge elements for a and 0-order surface elements for ϕ red . To describe the behavior of the magnetic field at the boundary, we add to this system of equations a discretized form of Green's third identity applied to the reduced magnetic scalar potential ϕ red . Indeed, at the boundary of the free space domain Ω 0 , ∆ϕ red = 0. Green's third identity applied to ϕ red reads,

1 2 ϕ red + ∂Ωm G B • n µ 0 d∂Ω - ∂Ωm ϕ red ∂G ∂n d∂Ω = ∂Ωm G H 0 • n d∂Ω. (8)
with G the Green's kernel. Equation ( 8) is then projected onto surface test functions constant by surface cell element, ϕ red is discretized with the same 0th-order surface shape functions. Magnetic induction leakage, B n , is linked to the circulation of the magnetic vector potential on the edge elements of the active domain boundary ∂Ω m , via an incidence matrix resulting from the local application of Stokes' theorem to the edges of a facet element. The outside and inside magnetic problems are then coupled via the continuity of B •n and ϕ red across ∂Ω m . BEM matrices are full and are compressed using the Fast Multipole Method (FMM) [START_REF] Greengard | A fast algorithm for particle simulations[END_REF] to reduce storage space up to 80

IV. SOLVING THE DISCRETIZED ME PROBLEM

The presence of both full BEM matrices and sparse FEM matrices makes the resolution of the fully coupled discrete problem difficult. This is associated with poor conditioning of the global system which results from the difference in the scale of the physics involved (MPa vs Fm -1 vs mH -1 ). A coupled multiphysics solution can be obtained from the iterative resolution of single-physics problems [START_REF] Urdaneta-Calzadilla | A FEM-BEM coupling strategy for the modeling of magnetoelectric effects in composite structures[END_REF].

An adaptive algorithm for the choice of the problem to be solved, which solves in priority the problem with the worst convergence, is used. MUMPS direct solver is used to solve the discretized electric equation while the Newton-Raphson method is used to solve the nonlinear mechanical and magnetic problems. They are based on a preconditioned BiCGSTAB solver for the discretized mechanical equation and a preconditioned GMRES for the discretized magnetic problem. An incomplete LU is used for the FEM block of the magnetic problem, while a Jacobi preconditioner is used for the diagonal BEM block.

V. VALIDATION OF THE MAGNETO-MECHANICAL

FORMULATION

The magneto-mechanical formulation was validated with respect to the analytical solution of a Galfenol sphere under a uniform source field and uniaxial stress. It involved solving for (B, S) in H = 3 H0 µr+2 , with µ r = B/µ 0 H, both supposed scalar and uniform because of the considered spherical geometry, with H 0 oriented along the (z) direction, together with solving for (B, S) in T zz = T 0 . Figure 1 shows the L 2 error of the formulation compared to the analytical solution [START_REF] Urdaneta-Calzadilla | A FEM-BEM coupling strategy for the modeling of magnetoelectric effects in composite structures[END_REF]. They show the h-convergence of the formulation for H 0 = 200 kAm -1 and T 0 = 10 MPa.

VI. ROTATING COILLESS ME DEVICE

The directions of easy magnetization in magnetostrictive materials are dependent on both the magnetic state and the state of stress. Depending on how these easy magnetization directions are controled in the magnetostrictive layer, two modes of operation of the ME composite are possible.

To illustrate these behaviors, we consider the structure presented in Figure 2 which constitutes a rotating coilless ME device. The device consists of a layer of a piezoelectric material, PZT-5A [START_REF] Jemai | Advanced Parametric Analysis of piezoelectric actuators with Interdigitated Electrodes having Various Cross-Sections[END_REF], on which four electrodes are arranged, and a layer of a disk-shaped magnetostrictive material consisting of Galfenol. This composite is subjected to the magnetic field of a NdFeB permanent magnet, with a remanent field of 1.13 T, free to rotate and spatially referenced by the angle θ (Fig. 2) corresponding to the direction of its polarization. An Amperian description is considered for the permanent magnet. The poling pattern of the piezoelectric phase is considered fixed, and the piezoelectric material poled by the application of electric potentials to a specific configuration of the electrodes: 1 kV on e 3 and e 4 and 0 V on e 1 and e 2 . This results in a poling pattern mainly along the θ = π/4 (Fig. 2). Numerically, this poling pattern was obtained by the rotation of the piezoelectric tensors C, e and ε in the direction of the electric field obtained by an electrostatic resolution.

A. Discrete system

The geometry was discretized into 27,055 nodes, 54,944 tetrahedral elements and 33,992 triangular elements for the surface of the magnet. This translates into 79,573 magnetic DoFs, 34,443 mechanical DoFs and 11,011 electrical DoFs. For every resolution of the mechanical problem, it takes around 21 iterations of the mechanical Newton-Raphson for the residue to converge within a relative tolerance of 10 -7 . The magnetic Newton-Raphson converges in 2 iterations to the same tolerance. The number of resolutions needed to reach a convergence of all solutions to relative 10 -6 in generator and actuator mode are presented in Table I and the convergence of the three single-physics solutions presented in Figure 3.

The total integration time of BEM and linear FEM matrices was of 36 s and the resolution time of each simulation of ME device took around 2 h on a computer equipped with Intel© Xeon© Gold 6240R CPU @ 2.4 GHz, 256 Gb of RAM.

B. Actuator Mode

In actuator mode, it is the stress field in the magnetostrictive layer, generated by the piezoelectric layer, which controls the directional magnetic anisotropy. By applying a positive electric voltage V in between the electrodes e 1 -e 2 and electrodes e 3 -e 4 , the piezoelectric layer mainly produces tension in the Galfenol layer. This induces a change in the direction of magnetization. As a result, the magnetic interaction between the Galfenol and the magnet produces a volumic torque, tending to align the polarization of the magnet with the easy axis of the Galfenol disk.

The torque Γ was computed by two methods. In the Galfenol layer (Ω mm ), the classical magnetic charge method [START_REF] Durand | Magnétostatique[END_REF] as been used:

Γ = Ωmm r × (∇ • M ) H magnet dΩ, (9)
with M the magnetization of the Galfenol layer and H magnet the magnetic source field created by the magnet. For the computation of the torque acting on the permanent magnet, a most efficient way is to consider equivalent current for the magnet magnetization. The torque is then simply computing considering the Laplace field current interaction.

Fig. 4. Global torque Γ between the Galfenol disk and the magnet vs the voltage applied to the piezoelectric layer. It is computed as a Laplace force in the magnet, and, its (z) component is compared to the torque computed by the magnetic charge method in the Galfenol layer.

Figure 4 shows the global torque as a function of the input voltage. To study the dependence of the angle θ onto the torque, the rotation of the magnet on a complete turn is discretized into 32 angles, and the magnetic force computed with the two previous methods for an input voltage V in = 1 kV. Each computation of the torque by the magnetic charge method takes around 22 s while the computation of the Laplace force takes around 3 h, for 64 Gauss point at each of the magnets triangular surface elements.

In Figure 5 we see the mechanical torque as a function of the angle θ. We see that for angles π/4, 3π/4, 5π/4 and 7π/4 rad, the torque is equal to zero. This is because at these angles the mechanical anisotropy induced by the mechanical loading of the Galfenol phase is in the same direction as the magnetization. There is therefore no global shifting in the direction of magnetization in the Galfenol phase and, therefore, the torque in the (z) direction is equal to zero. Conversely, angles 0, π/2, π and 3π/2 rad correspond to maxima of the torque. Given the poling pattern of the electrodes, which results in a poling and mechanical traction along the π/4 rad direction, the anisotropy of the magnetization in the Galfenol layer is maximal for θ equal to directions multiple of π/2 rad, i.e, π/4 rad relative to the π/4 rad main poling direction.

C. Generator mode

In the generator mode, it is the rotation of the magnet that controls the direction of the magnetization in the magnetostrictive layer and, thus, the deformations of the ME composite. This leads to the appearance of a electric potential difference between the electrodes, which is dependent on the angle of rotation θ (Fig. 2). The rotation of the magnet during a complet revolution is again discretized into 32 angles, and the previous poling pattern is considered. The electrodes are configured as follows: reference potential for electrodes e 1 and e 2 , floating potential for electrodes e 3 and e 4 . 6d show the amplified displacements (×10 5 ) of the device. They illustrate its mechanical response to the magnetic excitation. They show that the deformation of the device occurs mainly along the angle θ even if some slight anisotropies in the mechanical response are visible for θ = 0 & π/2 rad. These anisotropies in mechanical response can be explained by the anisotropic mechanical properties poled PZT-5A.

In Figure 7 is presented the output voltage of the ME composite vs the angle θ between the magnet and the structure. We observe a sinusoidal output voltage with two periods per rotation of the magnet. This is because the magnetostrictive strains depend only on the direction of the applied field and not on its sense.

VII. CONCLUSION

In this paper we present a FEM-BEM coupling strategy used for modeling nonlinear ME effects in composite structures. While using a formulation similar to [START_REF] Urdaneta-Calzadilla | Modeling of Magnetoelectric Effects in Composite Structures by FEM-BEM Coupling[END_REF], nonlinear coupled constitutive laws are used for the magnetostrictive phase. We adopted an invariant description of the Helmholtz free energy to derive the nonlinear expressions of the Cauchy stress tensor and the magnetic field in order to model the magnetomechanical behavior of Galfenol. We solved for the complete multiphysics problem by solving iteratively the three derived single-physics problems. The full algorithm was successfully tested on the modeling of a rotating coilless ME device in generator and actuator mode. In particular, two methods for computing the mechanical torque were presented and gave very close results. However, the magnetic charge method proved to be much faster at computing magnetic forces.

Fig. 1 .

 1 Fig. 1. L 2 error of the magneto-mechanical solutions compared to the analytical solutions vs the number of total DoFs of the magneto-mechanical problem.

Fig. 2 .

 2 Fig. 2. Test device, the arrow corresponds to the direction of the remanent field in the magnet, e 1 , e 2 , e 3 and e 4 refer to electrodes. Coefficients for Galfenol at 0 MPa has been obtained by a numerical fitting of experimental data from [5] : α 1 = 1514.2433, α 3 = 91.2778, α 5 = 64.8796, α 6 = -113.9582, α 7 = 44.2302, α 10 = -1.3699, α 11 = 0.2859, β 1 = -6.9705 • 10 6 , γ = 1.1280 • 10 10 , γ 2 = -6.9770 • 10 5 , α 2 = α 4 = α 8 = α 9 = 0.

Fig. 3 .

 3 Fig. 3. Convergence of the single-physics solutions in generator mode within the multi-physics solver for the ME generator at θ = 0 rad.

Fig. 5 .

 5 Fig. 5. Torque Γ vs rotation angle θ. Cross markers correspond to the z-component of the torque in the Galfenol layer (-Γ, z), circle markers correspond to the z-component of the torque in the magnet, for 33 angles θ describing the full rotation of the magnet.

Fig. 6 .

 6 Fig. 6. Amplified displacements (×10 5 ) of the ME composite for angles θ = 0, π/4, π/2 & 3π/4 rad, top view Figures 6a-6d show the amplified displacements (×10 5 ) of the device. They illustrate its mechanical response to the magnetic excitation. They show that the deformation of the device occurs mainly along the angle θ even if some slight anisotropies in the mechanical response are visible for θ =

Fig. 7 .

 7 Fig. 7. Potential difference between the reference electrodes, e 1 and e 2 , and the floating potential electrodes, e 3 and e 4 , vs the angle θ between the magnet and the ME composite.

TABLE I NUMBER

 I OF SINGLE-PHYSICS RESOLUTIONS NEEDED FOR CONVERGENCE OF THE MULTIPHYSICS PROBLEM WITH TOLERANCE INFERIOR TO 10 -6

	Problem	Generator Actuator (1 kV)
	Magnetic	4	4
	Electric	16	15
	Mechanical	16	15