B El Bahi

M Galai

M Cherkaoui

H Takenouti

Electrochemical deposition mechanism of copper-zinc-tin alloy and structural characterization

Keywords: Cu-Sn-Zn, Coating, Electrodeposition, Kinetics, Voltammetry, SEM, XRD

The electrodeposition of Cu-Sn-Zn coating (CZT) was investigated. The bath consisted of aqueous solution of Cu, Sn and Zn sulphate in the presence of sodium citrate (NaC6H5Na3O7) as complexing agent. The pH was set between 4 and 5.5 to stabilize the bath. Cyclic voltammetry exhibited four peaks during cathodic scan corresponding to the reductions of metal ions and adsorbed hydrogen. Impedance measurements were made at each potential peak. The effect of the metal salts concentration and the potential scan rate were studied. The deposit mechanism is controlled by mass transfer in parallel with a faradaic reaction. For the deposition of Cu-Sn and Cu-Sn-Zn alloys, the average diffusion coefficients of Cu 2+ and Sn 2+ ions and Cu 2+ , Sn 2+ and Zn 2+ ions were calculated. At E = -1.5 V, a Cu-Sn-Zn deposit is obtained with 55 at% copper, 25 at% tin and 20 at% zinc. The X-ray diffraction pattern shows the characteristic peaks of Cu, Cu6Sn5 and Cu3Zn2 phases. The deposit morphology was characterized by SEM.

Introduction

Thin-film solar cells have attracted the attention of numerous researchers in the last ten years. Initially, interest was focused on CuIn(Ga)Se2 (CIGS) due to their high energy conversion efficiency (21%) [1][2]. However, their development was slackened by the toxicity and the rarity of certain elements, namely indium and gallium [3][4]. Thus, kesterite structure alternatives, Cu2ZnSn (SxSe1-x)4 were proposed: Cu2ZnSnS4 (CZTS) [4][5][6][7][START_REF] Araki | Preparation of Cu2ZnSnS4[END_REF][9][10][11][12][13][14][15][16], Cu2ZnSnSe4 (CZTSe) [3,[17][18][19][20][21] and Cu2ZnSn(S,Se)4 (CZTSSe) [22][23][24][25][26][27]. These devices have high optical absorption coefficient (10 4 cm -1 ) and band-gap ranging between 0.9 and 1.5 eV, attractive for photovoltaic applications. Synthesis of CZTSe(S) layers, often used in these days, involves two steps: the formation of metallic precursor (CZT) and then the subsequent annealing under Se or S atmosphere at temperature between 450 and 600 °C. Different methods were used to develop the CZT metal precursor: spincoating [11,27], co-sputtering [24], pulse laser deposition [28], electrodeposition [21,25,[29][30], thermal evaporation [31][32], electroless plating [33], sol-gel method [34] and combustion [4]. The combustion method, based on the theory of propellant, is an extremely powerful tool to obtain metallic composite using corresponding metal salts-fuel mixtures. In some cases, it is necessary to combine two methods: galvanostatic deposition to fabricate a composition-spread film of binary Sn-Zn alloys, and then a second step involving a water gun with weak copper sulphate solution to achieve a thin film of ternary Cu-Sn-Zn alloys [35].

In this study, the formation of CZT film by electrochemical process is reported because this method is likely the most cost efficient, easiest to apply, and environmentally friendly. We noted that very few studies are devoted to the reduction mechanism of CZT alloy [26].

The simultaneous presence of Cu(II) and Sn(II) may make the electrolyte unstable because of a disproportionation reaction due to a high reducing power of stannous ions:

sol s sol sol Sn(IV) Cu Sn(II) Cu(II)    (1)
In order to obtain CZT alloy by electrodeposition, it is therefore necessary to add complexation agents in metal salts bath. Different agents were reported in the literature: citrate [10,12,26], tartrate [26],

or pyrophosphate [36]. In this study, we opted for citrate because it forms various electroactive aqueous complexes with Cu(II), Sn(II) and Zn(II) [26]. The starting bath composition was determined fromour previous work on the electrodeposition of nano-structured Sn-Zn alloy [37]. The concentration of copper sulphate and pH were fixed on the basis of the thermodynamic model proposed by Slupska and Ozga [38].

Fig. 1 presents the phase diagrams of dominant species of copper, tin and zinc in the presence of citrate (Cit) as complexation agent with respect to potential and pH. This figure was extracted from Fig. 4 of above quoted article [38]. In this figure, it can be remarked that at very cathodic potential (light green area), Cu-Sn-Zn alloy will be deposited. In the zone in cyan, Cu-Sn alloy is expected to form whereas in the brown area the metallic copper will deposit on the electrode surface. SnCit 2-will be stable in the white zone, but in this range Cu deposition will also take place. In the yellow zones (E ≥ ca. 0 V vs. SCE), hydroxide or hydroxyl-sulphate will precipitate rapidly, thus the bath composition in this area is not stable. Finally, in the light yellow zone, just above the stable SnCit 2-zone (more acidic medium), the solution is not entirely stable and precipitates may form slowly. Though not illustrated in this figure, for electrolytes, having a pH of more than 6.5, Sn(OH)2 will precipitate, thus the bath is not stable, neither.

Consequently, the electrode should be polarized in the zones where the metallic species, Cu or Cu containing alloys are stable [38], typically a solution with a pH between 3.5 and 6.5 and a potential less than 0 V vs. SCE.

Experimental conditions

.Materials and Preparation

The electrolysis cell was a borosilicate glass cylinder closed by a cap with five apertures. Three of them were used for the electrodes. The other two allow purging of the dissolved oxygen by nitrogen bubbling and temperature control. A copper with surface area of 1cm

Electrochemical methods

The electrochemical measurements were carried out using Voltalab PGZ 100® (Potentiostat/Galvanostat) monitored by Voltamaster 4.

Cyclic voltammetry

The voltammetric measurements were performed in the potential range between 0.6V and-1.8 V with a scan rate of 25 mV.s -1 .

Electrochemical impedance spectroscopy (EIS)

The impedance measurements were performed in a frequency range of 100 kHz to 10 mHz with 10 points per decade and with the AC amplitude of 10 mVrms. The EIS diagrams were presented as Nyquist plots. The results were then analyzed in terms of an equivalent electrical circuit.

XRD analysis

XRD spectra were collected with PANalytical X'Pert3 ® powder diffractometer with Cu Kα1 radiation. The semi quantitative deposit composition was determined by the Rietveld method using X'Pert high score ® plus software and CIF files from the American Mineralogist Crystal Structure Data base. The coatings morphology was examined with a FEI Quanta 200 ® scanning electron microscope.

Scanning Electron Microscopy (SEM)

The SEM images were carried using a (JOEL JSM-5500). The energy of the acceleration beam employed was 20 kV. Then, a steep current increase in absolute value occurs at -1.15 V due to the hydrogen evolution reaction.

Results and discussion

2H + + 2 e - H2 (3)
The citrate is thus electroactive species which will be reduced in the potential domain where the metal electrodeposition will take place. In presence of metallic cations (e.g. electrolyte 2), in the sodium citrate solution, four peaks can be observed during cathodic scan as illustrated in Fig. 3. This voltammogram is reported by various authors [26,[38][39].

XRD analyses

The species formed for each current peak were identified by XRD analyses.

The first peak at -0.30 V (a1) corresponds to Cu (II) reduction and is in agreement with the E-pH diagram of Fig. 1. Indeed, as can be seen in Fig. 4, the XRD pattern of the electrode surface confirms the presence of Cu in face centred cubic structure on the electrode surface [40][41][42]. Then, tin deposition starts around -0.70V and this wave reaches the maximum at -0.88 V (Fig 3, peak b1). Fig. 5 shows the diffraction pattern of the coating formed at this potential. We note the presence of cubic Cu, tetragonal tin and Cu6-Sn5 phases [18,34]. The third wave (c1) is attributed to reduction of adsorbed hydrogen (-1.1 V) and the fourth one (d1) to co-deposition of zinc together with Cu and Sn (-1.5V) confirmed by the presence of Cu3Zn2 phase diffraction peaks (Fig. 6) [43][44]. The deposition processes for Cu, Sn and Zn are in agreement with Fig. 1 [38] and also with the classification by Brenner [START_REF] Brenner | Electrodeposition of Alloys[END_REF].

Morphology and elemental distribution

Fig. 7 shows SEM images of deposits obtained at different potentials. We observe a granular, adherent and uniform coating at E = -0.3V (fig. 7a). It consists essentially of copper (Fig. 8). The granular structure is indicative of electrocrystallization under diffusion control conditions. At E= -0.88V, the deposits present some pores and irregularities due to the tin incorporation (Fig. 7b). The ternary coating Cu-Sn-Zn carried at -1.15 V is much more irregular and rougher (Fig. 7c). The EDX analysis confirmed that the composition of the deposit is related to the imposed potential. At E = -0.3 V, it is practically covered with copper alone whereas at E = -1.5 V, a Cu-Sn-Zn ternary deposit is obtained with 55at% copper, 25at% tin and 20at% zinc (Fig. 8).

Effect of bath compositions on the deposition kinetics

Effect of citrate concentration

Let us examine now the effect of bath compositions. Fig. 9 presents the effect of complexation agent concentration (Electrolytes 2, 3 and 4). In all cases, there are three current peaks for both cathodic and anodic potential scan, but the d peak at the highest citrate concentration was not observed because this peak is likely located outside the potential scanned domain. Fig. 10a shows that the potential of these peaks shift towards more negative values when the citrate concentration increases. It can be remarked also that the peak b is the most sensitive to the citrate concentration. The peak current increases as well as can be seen in Fig. 10b, when the citrate concentration increases. It should be recalled that the citrate reduction will take place at significant rate for a potential more negative than -0.7 V, i.e. in potential domain where peaks (a) to (d) were observed. The complexation of a metallic ion by a ligand, L occurs according to:

M n+ + L (ML) n+ (4)
With the complexation constant Kc, one gets:

(5)

The deposition of the metal M occurs according to:

(ML) n+ + ne- M + L (6)
With an equilibrium potential given by Nernst equation:

      L ML ln 0       n F n T R E E (7)
where E 0 is the standard potential of the (ML) n+ /M couple, n the number of exchanged electrons in the deposition reaction, F the Faraday constant (96487 C. mol -1 ), T the temperature (K), R the gas constant (8.31 J mol -1 K -1 ), [M n+ ] metal ion concentration, [L] ligand concentration and [(ML) n+ ] is complex concentration When the ligand concentration increases (citrate), the equilibrium (4) shifts to the forward direction leading to a decrease of M n+ concentration. Consequently, the equilibrium potential, according to Eq. 7 will decrease. Furthermore, we noticed an increase in the cathodic current density with the citrate concentration.

Effect of pH

As predicted by thermodynamic analyses by Slupska and Ozga [38], the bath is not stable when the bath pH was adjusted to 3. In contrast, it is stable between pH 4 and 5.5. Though not illustrated here, the shape of voltammograms obtained in this pH range is nearly unchanged. Fig. 11 summarizes the peak potential (Fig 11a) and the peak current density (Fig. 11b) with respect to the bath pH. The peak potentials, a, b and d decrease slightly when the solution pH becomes greater whereas the peak current densities increase when the bath becomes more alkaline.

The predominant forms of metal-citrate complexes are Cu2HCit2 3-, SnCit 2-, ZnHCit - [38], therefore the reactions taking place will likely be:

Cu 2 HCit 2 3-+ 3H + + 4e -2Cu + 2H 2 Cit 2-(8) SnCit 2-+ H + SnHCit - (9a) SnHCit -+ H + + 2e- Sn + H 2 Cit 2- (9b) ZnHCit -+ H + + 2e -Zn+H 2 Cit 2- ( 10
)
With equilibriums potentials for these three peaks:

(11)

(12) (13)
An increase in pH, i.e. a decrease in H + concentration induces the shift of the equilibriums processes (8 to 10) towards the left, then Eqs 11-13 confirm that pH increase shifts the peak potentials towards more negative direction. Eqs 8-10 may indicate that the increase of H + concentration increases the cathodic current, in contradiction with the experimental results. More detailed study is therefore necessary. Nonetheless, such behavior can be seen for pH 5 and 6.

Effect of the metal salts concentration

Figs 12 to 14 show the peak potentials and the peak current densities determined from the voltammograms at different concentrations of copper, tin and zinc sulphates respectively. The corresponding electrolyte is indicated in each caption. For these three metal ions, the concentration has little effect of the peak potential Ep. Nevertheless, the increase of metal salt concentration induces an increase of the cathodic current density.

Influence of the potential scan rate

To get information about the electrode kinetics, the voltammograms were obtained at different potential scan rates as displayed in Fig. 15.The X-axis is here the square root of the scan rate, v 0.5 . As can be seen in Fig. 15a, the peak potentials shift linearly with respect to v 0.5 towards more negative direction when the scan rate increases.

The peak current (Fig. 15b) varies also linearly. For the peak a, the line crosses the origin, therefore, the reaction rate is entirely controlled by the diffusion [37,[START_REF] Bard | Electrochemical methods; Fundamentals and Applications[END_REF]. In contrast, for the peaks b and d, the extrapolation to v = 0 give bias currents equal to -9 and -14 mA cm -2 respectively. The bias current at the peak d is probably due to the hydrogen evolution process as can be seen in Fig. 2, in parallel with the Cu and Zn reductions. The nature of the bias current of the peak b is the reduction of Cu controlled by the concentration gradient. Indeed, we deposit at E = -0.88V the Cu-Sn alloy and at E = -1.5V the Cu-Sn-Zn alloy.

The peak current density Ip corrected for the bias current for semi-infinite diffusion can be expressed as [START_REF] Bard | Electrochemical methods; Fundamentals and Applications[END_REF]: Where n represents the number of electron exchanged, D, diffusion coefficient (cm 2 s -1 ) v, potential scan rate (V s -1 ), and C, concentration of reacting species (mol cm -3 ).

  2 / 1 2 / 1 3 / 2 5 10 69 . 2 v C D n I p       (14)
Diffusion coefficient of copper (DCu 2+ ) is calculated from the line "a" of figure 15b. DCu 2+ = 1.45 10 -6 cm 2 s -1 , with ncu= 2, For the deposition of Cu-Sn and Cu-Sn-Zn alloys, it would be more appropriate to take into account the average diffusion coefficients of Cu 2+ and Sn 2+ ions and Cu 2+ , Sn 2+ and Zn 2+ ions calculated from the line "b" and "c" of figure 15brespectively. DCuSn = 1.29 10 -6 cm 2 s -1 , with ncu= nsn= 2, DCuSnZn = 1.41 10 -6 cm 2 s -1 , with ncu= nsn= nZn= 2 These values of diffusion coefficients are in agreement with some data reported in the literature when the deposition baths contain citrate [26,[START_REF] Ying | [END_REF]. In aqueous solution of infinite dilution, a value of 7.1 10 -6 cm 2 s -1 is reported for Cu and 3.8 10 -6 cm 2 s -1 for Sn [48][START_REF] Vanisek | CRC Handbook of chemistry and physics 91[END_REF].

Electrochemical Impedance Spectroscopy (EIS)

In Fig. 16, the red curve (-○-) shows the experimental impedance spectrum collected at -0.30 V corresponding to the reduction of Cu(II). The impedance spectrum exhibits a loop at high frequency attributed to the double layer capacitance in parallel with the charge transfer resistance Rt followed by the Warburg impedance in agreement with Fig. 15b. Indeed, the latter showed clearly that the electrochemical process is controlled by the mass transfer.

  Cd a Cd d t s Q Z R R Z        j 1 1
where

  W a d C D n Z         j 1 F T R 2 2 (15)
Where Rs, Rt, Zd, QCd, aCd, aW, and  represent respectively solution resistance ( cm²), charge transfer resistance ( cm²), diffusion impedance ( cm²), CPE parameter (F cm -2 s -1 ), CPE coefficient, and coefficient allowing a smaller slope of the Warburg impedance Zd ( cm²). The divergence from the unit slope for the Warburg impedance may be explained by a rough surface of electrodeposited copper. The result of parameter regression calculation carried out with house-made software (Simad, LISE) is presented with bleu curve (-+-). The agreement between the two curves is fairly good.

In Fig. 17, the red curve (-○-) is the impedance spectrum obtained at -0.88 V where the co-deposition of Cu and Sn is taking place. Below 0.1 Hz, a marked decrease of the impedance modulus can be noticed. This unusual feature was checked by a Kramers-Kronig transform. This transformation was carried out to determine the imaginary part from the real part in the admittance [START_REF] Gabrielli | Kramers-Kronig transformation in relation to the interface regulating device in « Electrochemical impedance: Analysis and interpretation[END_REF][START_REF] Barsoukov | Impedance Spectroscopy; Theory, Experiment, and Applications[END_REF]. The spectrum obtained from the Kramers-Kronig transform is presented by blue curve (-+-). A significant divergence was observed for the frequency below 1 Hz. A sharp decrease of impedance modulus is thus not validated and is due to the fact that the impedance spectrum was collected before reaching the steady-state conditions. The green curve (-x-) in Fig. 17 illustrates the impedance spectrum obtained by a non-linear regression calculation (between 1 kHz and 0.1 Hz) according to Eq. 14. No diffusion impedance can be seen clearly in the experimental spectrum (red curve), but Fig. 15b shows that the reaction rate is controlled by mass transfer with a bias faradaic reaction taking place in parallel.

The red curve (-○-) in Fig. 18 shows the experimental impedance spectrums collected at -1.50 V. As this was the case above, a marked decrease of the impedance modulus is noticed below 0.16 Hz. This feature was checked by the Kramers-Kronig transform as can be seen with blue curve with (-+-) symbol. The decrease of the impedance modulus is not validated, and then the impedance obtained was fitted with Eq. 15. The fitting was carried out between 3.16 kHz and 0.316 Hz, i.e. two extremities of the impedance spectrum were discarded.

The result of regression calculation is displayed by green curve with (-x-) symbol. In spite of poor experimental data, the fitted curve reproduces rather well the experimental spectrum between 1 kHz and 0.1 Hz.

R R R R Q C              1 1
The value of weighted with the impedance modulus both for the real and the imaginary part of the impedance, is rather small attesting suitable fitting quality.

The solution resistance Rs decreases when the applied DC potential becomes more cathodic. This variation is rather unusual and its origin is unclear. More detailed investigation in the electrolyte composition near the electrode surface will be necessary.

One of the tentative explanations is as follows: the water reduction taking place at the electrode surface makes the solution pH at the vicinity of the electrode surface smaller, i.e. [H + ] higher, thus the electrolyte conductivity may increase [53].

The charge transfer resistance Rt decreases when the DC current increases in the absolute value, except for the peak d. Cd is rather high, around 1 mF cm -2 , which can be explained by granular rough surface of electrodeposited alloys as can be seen in Fig. 7.

The value of aW is close to unity, that is, in spite of rough surface feature of the electrode, the distance between the top and the bottom of the surface might be small compared with the thickness of the diffusion layer. Note that the Warburg impedance was applied for fitting procedure, that is, the thickness of the diffusion layer is postulated to be much greater than that of ac concentration profile in the frequency range used.

The calculated concentration C  is also small compared with the concentration of metal cations added in the electrolyte. This might be explained by the complexation of diffusion species by citrate.

Conclusion

Cu-Sn-Zn coatings from citrate bath were successfully obtained. At -1.5 V vs. SCE, the deposited is composed of 55 at% of copper, 25 at% of tin and 20 at% of zinc. The X-ray diffraction pattern shows the characteristic peaks of copper, Cu6Sn5 and Cu3Zn2 phases. The voltammetric study showed that the deposition kinetics of Cu-Sn-Zn alloy is not controlled by pure diffusion unlike that of copper. Electrochemical impedance measurements confirm these results. For the copper deposit at E = -0.30 V vs. SCE, the impedance spectrum exhibits a loop at high frequency followed by the Warburg impedance showing that the process is controlled by the mass transfer. The experimental and fitted diagrams are in good agreement. At the deposition potential of the Cu-Sn-Zn alloy, the reaction rate is controlled by mass transfer with a bias faradaic reaction taking place in parallel. The fitted curve reproduces rather well the experimental spectrum between 1 kHz and 0.1 Hz.

Peak

Figure 1 :

 1 Figure 1: E-pH diagram of 0.01 M CuSO 4 , 0.1 M SnSO 4 , 0.25 M, ZnSO 4 , 0.65 M, Na 3 HCit [38]

3. 1 .Figure2:

 1 Fig. 2 displays the voltammogram obtained in electrolyte (1) containing sodium citrate only at pH = 4. The increase in the cathodic current density from -0.8 V is attributed to the reduction of protonated form of citrate ion according to the reaction: 2H3Cit -+ 2 e -H2 + 2H2Cit 2-(2)

Figure 3 :

 3 Figure 3: Cyclic voltamogram on a Cu electrode in 0.16 mol.L -1 CuSO 4 + 0.14 mol.L -1 SnSO 4 + 0.22 mol.L -1 ZnSO 4 and 0.4 mol.L -1 citrate at pH = 4 (electrolyte 2), scan rate v = 25 mVs -1

Figure 4 :

 4 Figure 4: XRD pattern of electrodeposited species in Electrolyte 2 at -0.30 V corresponding to Cu in cubic structure

Figure 6 :

 6 Figure 6: X-Ray diffraction pattern of electrodeposited coatings at -1.5V (electrolyte 2)

Figure 7 :

 7 Figure 7: SEM micrographs of deposits performed at : (a) E= -0.3V/SCE, (b) E = -0.88 V/SCE, (c) E= -1.5V/SCE.

Figure 8 :

 8 Figure 8: Deposit composition at three current peaks The peak position and the peak current density with respect to citrate concentration are illustrated respectively in Figs 10a and 10b. No correction for ohmic drop term was made.

Figure 9 :

 9 Figure 9: Cyclic voltammogramwith Cu substrate in electrolytes 2, 3, and 4, scan rate v = 25 mVs -1 .

Figure 10 :

 10 Figure 10: Effect of citrate concentration on (a) the peak potential E p, without ohmic-drop correction and (b) the peak current density I p.

Figure11:

  Figure11: Effect of bath pH on (a) the peak potential E p and (b) the peak current density I p (Electrolytes 2, 5, 6). Scan rate v = 25 mV s -1

Figure 12 :Figure 13 :Figure 14 :

 121314 Figure 12: Effect of CuSO 4 concentration on (a) the peak potential E p and (b) the peak current density I p (Electrolytes 2, 7, 8). Scan rate v = 25 mV s -1

Figure 15 :

 15 Figure 15: Effect of square root of scanning speed potential on (a) the peak potential E p and (b) the peak current density I p (electrolyte 2)

Figure 16 :

 16 Figure 16: Nyquist plot of EIS at -0.30 V (peak a) corresponding to the Cu depositionThis impedance spectrum was fitted with the following equation:

 cm 2 Cu+SnFigure 17 :

 217 Figure 17:Impedance spectrum collected at -0.85 V corresponding to CuSn deposition

 cm 2 Cu+Sn+ZnFigure 18 :

 218 Figure 18: Impedance spectrum collected at -1.50 V corresponding CuSnZn deposition

Table 1 :

 1 2 was used as working electrode; Pt plate as the counter electrode and a SCE as the reference electrode. Baths composition

	The commercial copper used in this study
	had the following chemical composition
	(by weight %): 0.019 P, <0.001 Fe, <0.001
	As, <0.001 Mn, <0.002Sb, <0.001 Al,
	0.009 Sn, 0.003 Ni, 0.015 Pb, <0.005 Ag,
	<0.001 Bi, <0.001 S, <0.005 C, the balance
	being Cu.
	The bath compositions examined are
	presented in Table 1, these solutions were
	blue, clear, and no precipitates were
	observed. The electrodeposition baths were
	made of copper sulphate (CuSO4•5H2O),
	tin sulphate (SnSO4), zinc sulphate
	(ZnSO4•7H2O), and trisodium citrate

.u. 2  / degree (c) E = -1.50 V vs. SCE

	5000
	2000 3000 4000 Intensity / aCu 6 Sn 5 (004) Cu 6 Sn 5 (012); SnO(012); Cu (111) CuSn(110); Cu 3 Zn 2 (033) SnO (132) SnO(111); Cu(0202) SnO(011) Cu 6 Sn 5 (006);
	1000	(00-3) Cu 6 Sn 5	Cu 3 Zn 2 (022)
	20 0	40	60	80	100

Table 2

 2

		summarizes the parameters
	obtained by regression calculations. In this
	calculation, the diffusion coefficient was
	fixed because the calculation cannot be
	separate D and C .
	In this table, instead of QCd and aCd, the
	effective capacitance according to Brug et
	al. was used [52].
						a	Cd
	d	Cd	Cd a	s	t	aC	d
				s	t

Table 2 :

 2 Kinetics parameters determined by the fitting calculation (Eq. 14).

Acknowledgements

The authors wish to thank Dr. E. Chassaing (Institute of Research and Development on Photovoltaic Energy, EDF-CNRS-Chimie Paristech ENSCP) for her pertinent remarks on this manuscript.


domain: spi phys math info

text:
Scott Delbecq

Marc Budinger

Clément Coic

Nathalie Bartoli

Nathalie Bartoli Trajectory

Trajectory and design optimization of multirotor drones with system simulation

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

The design optimization of multirotor drones is a key enabler for improving existing and future vehicle architectures performance. For many multirotor applications, the design is performed in order to respect the maximum speed requirements and maximizing the endurance in hover of the drone. For specific applications, the requirements enable to propose different alternatives for the trajectory which can thus be considered as an additional degree of freedom in the design process. In this paper, different approaches are proposed for achieving the vehicle design optimization, the trajectory optimization and both simultaneously. The approaches taken explore the possible use of multidisciplinary analysis and design optimization with Open-MDAO and dynamic optimization with Modelon Impact and Optimica. The specific multirotor drone application considered is a lifting and handling vehicle which could replace cranes in dense urban areas. The most significant finding is that vehicle operations that involve transient trajectories considerably affect the overall drone design especially the choice of the propeller pitch. = Oversizing coefficient = Mass = Multidisciplinary Design Optimization * Researcher, Department of Aerospace Vehicles Design and Control, scott.delbecq@isae-supaero.fr. † Associate Professor, Department of Mechanical Engineering, INSA Toulouse, marc.budinger@insa-toulouse.fr. ‡ Senior Simulation Engineer, Dr.-Ing., Industry Solution Services, clement.coic@modelon.com. § Senior Researcher, Information Processing and Systems Department, nathalie.bartoli@onera.fr, AIAA Member.

Nomenclature

I. Introduction

T design optimization of Unmanned Aerial Vehicles (UAV) has gained interest due to the variety of applications they offer. Multirotor drone applications cover both civil and military applications including surveillance, agriculture and other unprecedented purposes, such as automated package delivery or Personal Air Vehicles (PAVs) because of their high maneuverability and capability to hover. However, their endurance is their bottleneck when compared to winged architectures. If a wide range of applications is sought such as for consumer multirotor drones, the vehicle design is performed for a given maximum forward flight speed and such that the the hovering time or range in forward flight is maximized [START_REF] Kim | Development of a Comprehensive Analysis and Optimized Design Framework for the Multirotor UAV[END_REF]. The analysis [START_REF] Shi | A practical performance evaluation method for electric multicopters[END_REF] and design optimization [START_REF] Bershadsky | Electric multirotor UAV propulsion system sizing for performance prediction and design optimization[END_REF][START_REF] Delbecq | Efficient sizing and optimization of multirotor drones based on scaling laws and similarity models[END_REF] of multirotor drones have been thoroughly achieved for such requirements. For specific applications, the vehicle design may consider a specific mission for which the vehicle is tailored perfectly resulting in optimized mass and cost but limited usage. Often a mission can be achieved using different trajectory alternatives with regard to path and/or time. Thus, the trajectory may be considered as an additional degree of freedom in the design process for specific applications. Such approach has been used for winged vehicles design optimization [START_REF] Falck | Trajectory optimization of electric aircraft subject to subsystem thermal constraints[END_REF] but not multirotor drones and their applications. In addition, comparison of strategies such as distributed or monolithic for solving the trajectory and design optimization problem are absent in the literature. In this paper, we aim at achieving the trajectory and design optimization of a lifting and handling point-to-point [START_REF] Park | Motion profile planning of repetitive point-to-point control for maximum energy conversion efficiency under acceleration conditions[END_REF] multirotor drone.

Two computation environments are used for implementing the models and solving the optimization problems. The first one is FAST-OAD [START_REF] David | From FAST to FAST-OAD: An open source framework for rapid Overall Aircraft Design[END_REF] framework which is based on OpenMDAO [START_REF] Gray | OpenMDAO: An open-source framework for multidisciplinary design, analysis, and optimization[END_REF]. For this environment, the SLSQP optimizer is used and the semi-analytic finite difference approach is chosen for computing derivatives [START_REF] Delbecq | Benchmarking of monolithic MDO formulations and derivative computation techniques using OpenMDAO[END_REF]. The second one is Modelon Impact [START_REF] Coïc | Collaborative Development and Simulation of an Aircraft Hydraulic Actuator Model[END_REF] that uses the Modelica modelling language [START_REF] Elmqvist | Modelica-a language for physical system modeling, visualization and interaction[END_REF] and the Optimica language [START_REF] Åkesson | Optimica-an extension of modelica supporting dynamic optimization[END_REF][START_REF] Åkesson | Modeling and optimization with Optimica and JModelica. org-Languages and tools for solving large-scale dynamic optimization problems[END_REF]. The first advantage of Modelica language can deal with both Initial Value Problem (IVP) and Boundary Value Problem (BVP) for the model initialization. In the former case, the user provides the initial values for every state. In the latter, the number of independent initial values should match the number of states but is not necessarily their initial values. Thus, Modelon Impact compiler solves an initial problem -different from the dynamic problem -to resolve the BVP. In this process, the compiler can solve linear and non-linear systems, which is often not tolerated out of the box with a different solution. A second advantage is that the Modelica language is acausal. This way, defining a Drone model based on the equations of the physics makes it useful for several use cases: position, speed or torque command. In our example, we want to optimize the torque trajectory while sizing the drone. Nevertheless, it is convenient to provide the optimizer a start trajectory not too far from our constraints. This is easily achieved by simulating the same model, providing a position trajectory which not optimized at all but matching our requirements.

To model the trajectory a lumped parameter model is used and is implemented using the Modelica modeling language . The model is then exported from Modelon Impact using the Functional Mock-up Interface (FMI) standard [START_REF] Blochwitz | Functional mockup interface 2.0: The standard for tool independent exchange of simulation models[END_REF]. Optimica can directly control the FMU whereas to be controlled from Python and FAST-OAD the package pyFMI [START_REF] Andersson | Pyfmi: A python package for simulation of coupled dynamic models with the functional mock-up interface[END_REF] is used.

To study the trajectory and design optimization of a lifting and handling multirotor, the work is decomposed in three parts. First, the design optimization of the vehicle is achieved using a sizing model implemented in FAST-OAD using the Normalized Variable Hybrid (NVH) formulation [START_REF] Delbecq | Benchmarking of monolithic MDO formulations and derivative computation techniques using OpenMDAO[END_REF]. Then, the trajectory optimization problem is achieved on the previous design carried out. Lastly, the trajectory and vehicle design optimization is performed simultaneously using two approaches. First, the trajectory analysis is achieved by a Functional Mock-up Unit (FMU) and included in the vehicle overall design optimization in FAST-OAD. Finally, the whole problem is formulated in Modelica and Optimica using Modelon Impact.

To present the proposed work, the paper is organized as follows. In Section II, the lifting and handling multirotor case study is presented. Section III presents the vehicle design optimization problem and results. Section IV describes the trajectory model and optimization results. In Section V, the trajectory and vehicle design optimization is achieved and compared with the decoupled approach. Finally, Section VI offers concluding remarks.

II. Case Study: Lifting and Handling Multirotor Drone

Multirotor drones have become popular test-beds for UAV research and applications partly because of their great maneuverability and ability to execute complex movements. The multirotor capabilities have enabled them to be implemented in a wide range of applications. The components to be selected must therefore mainly allow movements and trajectories such as:

• Hover and low-speed movements, for video recording for cinema or security;

• Horizontal and vertical load transport for medical delivery or fire fighting;

• Surface scanning for agricultural spraying, Lidar Remote Sensing, area inspection;

• Complex missions assembling the different movements mentioned above. Different architectural choices can be made depending on the purpose and the mission. For multirotor drones, the main choices are the number of arms and the number of propeller per arms but also the materials and technologies of the components. In the presented work, a single architecture is considered for the multirotor drone. The architecture is presented in Figure 1 and is composed of:

1) 4 fixed pitch propellers 2) 4 out-runner brushless motors 3) 4 electronic speed controllers (ESC) mainly made from MOSFET inverters 4) 1 battery based on Li-Ion cells 5) A mechanical structure (frame) consisting of 4 arms and a central body Such vehicle could replace cranes for lifting at the top of a building some relatively small missing items or materials ( 25 ) when approaching the end of the construction. Some plausible reasons to replace the crane is to avoid a daily costly rent of a crane for a few items and also when the ground space remaining is small such as in urban areas as shown in Figure 2. In addition, this solution enables to avoid another solution which consists in employing interim workers to carry up the items with their natural strength and thus affecting their physical health.

This case study leads to interesting scientific challenges such as the design of an all-electric 25 payload multirotor and the determination of an optimal point-to-point trajectory.

The two main requirements for this case study are:

• Payload: = 25 . • Endurance: 150 climbs of height ℎ = 10 with 5 seconds of hovering for the handling of the payload. To show the interest in simultaneously determining the trajectory and the design, the vehicle is first design without usage of system simulation for the trajectory. Then, the analysis and optimization of the trajectory is achieved for the obtained design. The resulting trajectory and design will then be compared to the ones obtained with the holistic

III. Vehicle Design Optimization

A. Sizing Scenarios, Design Drivers and Models

The main sizing scenarios, design drivers and models for vertical flight applications of multirotor drones are summarized in Figure 3. Such applications consist of three sizing scenarios to be considered in the design problem that are the hovering flight with the advance ratio of the propeller = 0, the takeoff phase which requires maximum power to accelerate the drone with an increasing and the climb phase with a constant vertical speed and thus constant .

For the design optimization, the sizing scenarios considered are hover and takeoff which respectively contribute to the continuous power and maximum power sizing of the components. This approach enables to meet consistently endurance in hovering and takeoff acceleration capability objectives. The component sizing models used are scaling laws, linear regressions of data sheet and surrogate models. These are described in a thorough manner in [START_REF] Budinger | Scaling laws and similarity models for the preliminary design of multirotor drones[END_REF].

For the considered application, the lifting and handling system, the sizing scenarios have to be adapted. It is not straightforward to specify an optimal climb time or climb energy to be achieved. Hence for the design optimization, without system analysis for the trajectory analysis, it is proposed to specify a takeoff acceleration and height to determine the travel time using typical energy-efficient point-to-point mission profiles. Here, as the the inertial effects dominate the additional forces to those of the hovering flight, the triangular speed profile is chosen [START_REF] Reysset | Conception préliminaire d'actionneurs électromécaniques -outils d'aide à la spécification et à la génération de procédures de dimensionnement pour l'optimisation[END_REF][START_REF] Richiedei | Analytical computation of the energy-efficient optimal planning in rest-to-rest motion of constant inertia systems[END_REF]. The maximum acceleration for such profile is expressed as follows: where H is the displacement [START_REF] Coïc | Collaborative Development and Simulation of an Aircraft Hydraulic Actuator Model[END_REF] and the travel time. If we choose = 5 ( / 2 ) for takeoff , around the half of the gravitational constant , this yields = 2.8

To size the multirotor, the energy of the mission is obtained as follows:

= ,ℎ • ℎ + , • • (2)
Where is the number of displacements, ,ℎ and , are respectively the electric power during hover and takeoff at the battery level.

B. Problem Formulation

The overall sizing model has been adapted from [START_REF] Delbecq | Efficient sizing and optimization of multirotor drones based on scaling laws and similarity models[END_REF] which is tailored for vertical flight applications of multirotor drones. The formulation proposed in [START_REF] Delbecq | Efficient sizing and optimization of multirotor drones based on scaling laws and similarity models[END_REF] and illustrated in Figure 4 uses the Normalized Variable Hybrid formulation proposed in [START_REF] Delbecq | Benchmarking of monolithic MDO formulations and derivative computation techniques using OpenMDAO[END_REF]. This formulation is very efficient regarding computational cost and robustness to scale change and thus has been chosen for the lifting and handling system case study. The adaptation of the sizing scenarios leads to some changes in the problem formulation. The new design optimization problem can be formulated as follows: minimize 1 with respect to 1

β (0) pro , k (0) M , k (0) mot , k (0) speed,mot k (0) ND , k (0) mb , k (0) vb , J (0) , k (0) ESC ,
1 1 1 , 1 , 1 , 1 , 1 subject to , - , ≤ 0 1 - ≤ 0 1 - ≤ 0 1 - ≤ 0 1 - ≤ 0 1 ( 3
)
The objective is to minimize the weight of the vehicle for the defined mission with respect to design variables such as the propeller pitch or the propeller speed/diameter coefficient . The design has to respect some constraints to respect the technological constraints of components such as motor maximum electromagnetic torque , as well as others to respect voltages consistency in the power train ( ≤ and ≤ ). Some consistency constraints are used to solve multidisciplinary couplings as suggested by the NVH formulation [START_REF] Delbecq | Benchmarking of monolithic MDO formulations and derivative computation techniques using OpenMDAO[END_REF] ( ≤ and ≤

). This also requests to add normalized design variables such as and , .

C. Results

The problem has been solved using a Python script with the SLSQP algorithm of the Scipy package and using the Overall Aircraft Design framework FAST-OAD [START_REF] David | From FAST to FAST-OAD: An open source framework for rapid Overall Aircraft Design[END_REF] which uses OpenMDAO [START_REF] Gray | OpenMDAO: An open-source framework for multidisciplinary design, analysis, and optimization[END_REF]. The two tools converge to the same design with the main characteristics given in Table 1. The obtained design is rather tailored for hovering flight as the solution converges to a small propeller pitch . The propeller design converges towards relatively low rotational speeds with a large diameter. The vehicle mass is almost equal to the payload mass with more than one half of the mass coming from the battery. The design is capable of holding more than 28 minutes in hovering flight. The use of 2 propellers per arm in push-pull configuration could have allowed the use of higher rotational speed in order to gain motor mass.

IV. Trajectory Optimization

A. Models

For the multirotor drone design, it is important to be able to access the expressions of power variables at the different components levels. For the flight conditions in static hover, with no forward speed, no vertical speed, simple algebraic equations can be set up to assess the desired physical quantities. This scenario has a significant effect on the autonomy of the vehicle. The thrust generated by propellers have to compensate the overall weight of the drone and the payload. The propeller performances can be expressed thanks to aerodynamic coefficients and : the thrust = 2 2 , the power (and thus the torque) with = 3 5 with the propeller speed in rev/s, the diameter in m, the air density in kg/m3. It is then possible to calculate the flow of energy from the motor to the battery using the classical relationships of electrical components.

In the case of a usage scenario involving a more complex trajectory, such as the lifting and handling system,it is possible to use a model with lumped parameters like the one presented in Figure 5.

This model is implemented using the Modelica modeling language. It enables to analyze and solved the physics of a vertical flight consisting of the following Differential Algebraic Equations (DAE):

= (4)
= • (7) = • • 2 • 4 (8) = - 1 2 2 - (9) = • • 3 • 5 (10) = • 2 (11) ( + ) = - (12)
= (13) = + ( 14
) , = • ( 15
) , = , • (16)
= ∫ , (17)
This model enables to evaluate different variables of interests with respect to an input motor torque vector and different parameters obtained in the design process as shown in Figure 6.

B. Problem Formulation

Two tools were used to carry out the trajectory optimization. First, the Modelon Impact and the Optimica solver [START_REF] Åkesson | Optimica-an extension of modelica supporting dynamic optimization[END_REF][START_REF] Åkesson | Modeling and optimization with Optimica and JModelica. org-Languages and tools for solving large-scale dynamic optimization problems[END_REF] are used to solve the problem as an dynamic optimization problem. Then, the previous Modelica model is exported using Modelon Impact as a Functional Mock-up Unit (FMU) and integrated in OpenMDAO [START_REF] Gray | OpenMDAO: An open-source framework for multidisciplinary design, analysis, and optimization[END_REF] and FAST-OAD [START_REF] David | From FAST to FAST-OAD: An open source framework for rapid Overall Aircraft Design[END_REF].

The goal here is to evaluate the energy consumption ( ) to transport the payload load to a defined height so that we can assess the required battery autonomy for multiple elevations. The approach taken for the optimization trajectory problem is to use the motor torque vector of size in order to easily take into account the torque limits of the motors and the effects of the inertia of the propellers. The optimization problem can be formulated as follows: minimize 1 with respect to ( )

1 subject to ( ) - , ≤ 0 ℎ -( ) ≤ 0 1 ( ) = 0 1 (18)
The trajectory optimization problem consists in minimizing the energy for a single travel with respect to the motor torque time vector and the time to reach the desired position . Also, two mission profile constraints are used which are the height to be reached and a rest position ( ( ) = 0) in addition to a design constraint which is the motor peak torque. An initial and reasonable torque profile can be obtained by inverse simulation [START_REF] Delbecq | Knowledge-Based Multidisciplinary Sizing and Optimization of Embedded Mechatronic Systems-Application to Aerospace Electro-Mechanical Actuation Systems[END_REF][START_REF] Liscouët | Modelling approach for the simulation-based preliminary design of power transmissions[END_REF], which can simplify the work of optimization algorithms.

C. Effect of Travel Time, Maximum Torque and Number of Control points

The purpose of this part is to show the effect of different settings on the resulting trajectories. First, the effect of choosing the travel time as a design variable or fixed to 2.8 is studied. For this purpose, it is chosen to do the assessment for different amounts of maximum motor torque. The low torque value is close to the one found in the design optimization whereas the high torque is the double of this value. Four different optimization problems are thus performed which results in four different trajectories as illustrated in Figure 7.

The second study concerns the effect of the number on control points of the input torque vector ( ), hence the size of the design variable. Four different trajectory optimization are achieved with respectively = 4, = 8, = 16 and = 32. The resulting trajectories are given in Figure 8. As shown in Figure 9, increasing the number of control points enables to find a better optimum but increases the computational cost.

A number of control point greater than 16 does not lead to a significant decrease in the objective to minimize ( ). Therefore, for the trajectory optimization and the simultaneous trajectory and vehicle design optimization = 16 is chosen. In addition, the constraint on the trajectory for respecting the maximum torque of the motor has a significant effect on the resulting trajectory.

D. Results

Figure 10 illustrates the obtained trajectory for the design achieved in the previous section. The torque is limited by the motor peak torque constraints resulting in a moderate acceleration and thus longer travel time 3.45 . It can be depicted that the speed profile is not symmetrical as for triangular or sinusoidal speed profiles in order to decelerate at the latest moment and thus decrease travel time. Fig. 10 Resulting trajectory for the trajectory optimization approach using OpenMDAO and FMU.

V. Simultaneous Trajectory and Vehicle Design Optimization

A. Problem Formulation

The trajectory analysis is now added to the overall design optimization problem leading to an accurate estimation of the travel energy. Also, the trajectory optimization design variables and constraints are merged with the ones of the design optimization problem. The objective is the one of the design optimization problem that is the of the vehicle. A monolithic optimization approach is used and was solved both in the FAST-OAD and Modelon Impact environments. The problem is formulated as follow and consists of 10 design variables, 7 inequality and 1 equality constraint: minimize 1 with respect to 1

1 subject to , - , ≤ 0 1 - ≤ 0 1 - ≤ 0 1 - ≤ 0 1 - ≤ 0 1 ( ) - , ≤ 0 ℎ -( ) ≤ 0 1 ( ) = 0 1 1 1 1 , 1 , 1 , 1 , 1 ( )

B. Results and Discussions

The two environments converge to a similar design solution. However, the computation time is higher using FAST-OAD (> 2min) as the code is not compiled like with Optimica (< 1min). The computational cost in FAST-OAD may be potentially be reduced if the trajectory analysis is achieved using a surrogate model built with the Surrogate Modeling Toolbox (SMT) [START_REF] Bouhlel | A Python surrogate modeling framework with derivatives[END_REF] and analytic derivatives implemented when possible.

Figure 11 shows the resulting trajectory with the one obtained with trajectory optimization on the fixed design. One can see that the optimizer finds a better optimum by reducing the travel time to 2.89 at the cost of higher torque and mass motor. This also increases the mass of the other components and especially the battery as shown in Table 2. The pitch is increased for better vertical speed performance.

Nevertheless, the optimizer is able to find the best compromise between the travel time and intensity of the acceleration phases which indirectly correspond to a compromise between the weight of motors and the weight of the battery pack. One can observe that the initial prediction of the travel energy has been underestimated compared to the trajectory analysis. This results that the first design does not meet the 150 elevations requirements but has only the energy to achieve 100 of them. The MTOW of the vehicle is rather high but meets the initial requirements for this specific lifting and handling application. This work is a step forward for the integration of higher fidelity modules in FAST-OAD such as the sizing and vibration analysis of structural components using Finite Element Model simulations.

VI. Conclusions and Future Work

The work presented in this paper was built around the design and optimization of a lifting and handling multirotor. First, the specifications of such system were given. Then, a preliminary design was achieved based on an approximation of the mission energy. Then for the resulting design, a trajectory optimization was performed. Simulations showed the interest of setting the travel time as a degree of freedom and the importance of considering the maximum motor torque constraint defined by the design. Finally, the simultaneous trajectory and design optimization was performed.

The optimal solution tends towards reducing the travel time in the limit of what amount of torque the motor can provide. The trajectory and optimization results show that the dynamics of the trajectory affect the overall drone design especially the choice of the propeller pitch. Finally, predicting the optimal trajectory and specifying the energy to the design is not straightforward and underlines the need to achieved design optimization and trajectory optimization simultaneously as well as having an accurate trajectory analysis. The work also proved that this analysis can be achieved by integrating system simulation into MDAO process using FMUs.

Future work will consist in exploring different architectures for the multirotor such as push-pull configurations. Also, a detailed comparison and benchmark between the dynamic optimization approach of Optimica and MDAO approach of OpenMDAO for this type of design optimization problem could be of interest.

=

  Pitch/diameter ratio [-] = Propeller Thrust Coefficient [-] = Propeller Power Coefficient [-] = Diameter [ ] = Differential Algebraic Equation = Electronic Speed Controller = Energy [ ] = Force [ ] = Functional Mock-up Interface = Functional Mock-up Unit = Gravitational constant [ 3 -1 ˘2] ℎ = Building height [ ] = Advance ratio [-]

1 =

 1 Multidisciplinary Design Analysis and Optimization = MAximum Takeoff Weight = Rotational speed [ / ]

Fig. 1 Fig. 2

 12 Fig. 1 Multirotor Drone Architecture.

4 Fig. 3

 43 Fig. 3 Design drivers and sizing scenarios for vertical flight multirotor drone design [4].

Fig. 4 XDSM

 4 Fig. 4 XDSM diagram of the for vertical flight multirotor design optimization using the NVH formulation [4].

Fig. 5 Vertical

 5 Fig. 5 Vertical Flight Multirotor Drone Modelica Model.

Fig. 6

 6 Fig. 6 Parameters and variables of interest of the trajectory analysis.

Fig. 7

 7 Fig. 7 Resulting trajectories for fixed and free travel time with different torque amounts.

Fig. 8

 8 Fig. 8 Resulting trajectories for the different number of control points settings.

Fig. 9

 9 Fig. 9 Effect of the number of control points on the objective and the number of function calls.

Fig. 11

 11 Fig. 11 Comparison of resulting trajectories for the trajectory optimization approaches with and without the design.

Table 1 Resulting main characteristics for the design optimization of the multirotor.

 1

	Parameter	Value	Units	Description
		25.0	[kg]	Payload
		50.8	[kg]	Maximum Takeoff Weight
	ℎ ,	28.6	[min]	Endurance in hover
		0.5	[ 3 / / 2 ] Max equivalent acceleration ( ) during takeoff
		0.3	[m]	Propeller pitch
		1.134	[m]	Propeller diameter
	Ω ,ℎ	1660	[rpm]	Motor speed during hover
	,ℎ	7.1	[N.m]	Motor torque during hover
	Ω ,	2033	[rpm]	Motor speed during takeoff
	,	10.6	[N.m]	Motor torque during takeoff
		0.802	[m]	Arm length
		0.247	[kg]	Single propeller mass
		1.8	[kg]	Single motor mass
		15.4	[kg]	Battery mass
		1.8	[kg]	Frame mass

Table 2 Resulting main characteristics for the trajectory and design optimization of the multirotor.

 2

	Parameter	Value Traj. and Design Opt.	Value Traj. Opt.	Units	Description
		25.0	25.0	[kg]	Payload
		69.9	50.8	[kg]	Maximum Takeoff Weight
	ℎ ,	36.3	28.6	[min]	Endurance in hover
		150	100	[-]	Number of displacements to the top of the building possible
		2.89	3.45	[s]	Time to reach the top of the building
		0.61	0.5	[ 3 / / 2 ]	Max equivalent acceleration ( ) during takeoff
		0.43	0.3	[m]	Propeller pitch
		1.329	1.134	[m]	Propeller diameter
	Ω ,ℎ	1230	1660	[rpm]	Motor speed during hover
	,ℎ	17.3	13.1	[N.m]	Motor torque during hover
	Ω ,	1560	2033	[rpm]	Motor speed during takeoff
	,	21.1	10.6	[N.m]	Motor torque during takeoff
		0.939	0.802	[m]	Arm length
		0.339	0.247	[kg]	Single propeller mass
		3.2	1.8	[kg]	Single motor mass
		26.7	15.4	[kg]	Battery mass
		3.1	1.8	[kg]	Frame mass

VII. Acknowledgements

The authors would like to thank members of the GIS SizingLab for New Mobility and the ONERA -ENAC -ISAE-SUPAERO research federation.