Koussaila Hamiche
email: koussailahamiche@gmail.com

Michel Fliess
email: michel.fliess@polytechnique.edu

Cédric Join
email: cedric.join@univ-lorraine.fr

Hassane Abouaïssa
email: hassane.abouaissa@univ-artois.fr

Bullwhip effect attenuation in supply chain management via control-theoretic tools and short-term forecasts: A preliminary study with an application to perishable inventories

Keywords: Supply chain management, inventory control, perishable inventories, bullwhip effect, delay systems, Smith predictors, model-free control, short-term forecasting, time series, high gain

Supply chain management and inventory control provide most exciting examples of control systems with delays.

Here, Smith predictors, model-free control and new time series forecasting techniques are mixed in order to derive an efficient control synthesis. Perishable inventories are also taken into account. The most intriguing "bullwhip effect" is explained and attenuated, at least in some important situations. Numerous convincing computer simulations are presented and discussed.

I. INTRODUCTION

Supply chain management [START_REF] Bozath | Introduction to Operations and Supply Chain Management[END_REF], [START_REF] Dolgui | Supply Chain Engineering[END_REF], or SCM, i.e., the management of the flow of goods and services, is a salient subfield of operations research. Tools from control engineering have been tested since the very beginning (see, e.g., [START_REF] Simon | On the application of servomechanism the study of production control[END_REF], [START_REF] Vassian | Application of discrete variable servo theory to inventory control[END_REF]). Unavoidable dead times, often called lead times or throughput times, are due to phenomena like transport and production. Delay differential equations play therefore an increasing rôle (see, e.g., [START_REF] Abbou | On inventory control for perishable inventory systems subject to uncertainties on customer demands[END_REF], [START_REF] Escudero-Gomez | Inventory trajectory tracking via an demand-rate-estimation using Smith's predictor[END_REF], [START_REF] Chołdowicz | Comparison of a perpetual and PD inventory control system with Smith predictor and different shipping delays using bicritical optimization and SPEA2[END_REF], [START_REF] Dejonckheere | Measuring and avoiding the bullwhip effect: A control theoretic approach[END_REF], [START_REF] Hamiche | Decentralized control of dynamic supply chain systems with parametric uncertainties[END_REF], [START_REF] Ignaciuk | Smith predictor based control of continuous-review perishable inventory systems with multiple supply alternatives[END_REF], [START_REF] Ignaciuk | Dead-time compensation in continuousreview perishable inventory systems with a remote supply[END_REF], [START_REF] Moussaoui | Robust inventory control of production systems subject to uncertainties on demand and lead times[END_REF], [START_REF] Riddalls | The stability of supply chains[END_REF], [START_REF] Riddalls | Production-inventory system controller design and supply chain dynamics[END_REF], [START_REF] Schwartz | A process control approach to tactical inventory management in production-inventory systems[END_REF], [START_REF] Schwartz | A control-relevant approach to demand modeling for supply chain management[END_REF], [START_REF] Schwartz | Simulation-based optimization of process control policies for inventory management in supply chains[END_REF], [START_REF] Sipahi | On stability problems of supply networks constrained with transport delay[END_REF], [START_REF] Warburton | An exact analytical solution to the production inventory control problem[END_REF], and the references therein). Several previous publications use a model, that is more or less analogous to the following delay system, which is derived via a straightforward conservation principle,

ẏ(t) = ku(t -L) -d(t), t ≥ L (1)
where

• the output variable y(t) ≥ 0 is the material inventory,

• the control variable u(t) ≥ 0 is the factory supply,

• d(t) ≥ 0 is the customer demand,

• L ≥ 0 is the throughput time, which is known,

• k 0 is a yield parameter, which might be poorly known. The above infinite-dimensional system looks easy. The numerous corresponding control strategies in the literature do not seem however to be entirely satisfactory. Let us grasp the difficulty via the two following points:

• System (1) is marginally stable [START_REF] Franklin | Feedback Control of Dynamic Systems[END_REF].

• Assume, for instance, that we want to keep the inventory constant, i.e., ẏ(t) = 0. Equation (1) yields

u(t) = d(t+L) k
. The future value d(t + L) of the customer demand has to be "guessed." No theoretical technique will ever produce rigorously accurate predictions. For simplicity's sake let us presume here a constant bias b = 0, i.e., ẏ(t) = b k . Thus y(t) = b k (t -L) + y(L), t ≥ L. Its absolute value | y(t) | becomes very large, at least from a purely mathematical standpoint, when t → +∞. This is a convincing illustration of the famous bullwhip effect to which a huge literature has been devoted (see, e.g., [START_REF] Chen | Quantifying the bullwhip effect in a simple supply chain: The impact of forecasting, lead times, and information[END_REF], [START_REF] Dolgui | Supply Chain Engineering[END_REF], [START_REF] Lee | The bullwhip effect in supply chains[END_REF], [START_REF] Sucky | The bullwhip effect in supply chains -An overestimated problem?[END_REF], [START_REF] Udenio | Behavioral causes of the bullwhip effect: An analysis using linear control theory[END_REF], [START_REF] Wang | The bullwhip effect: Progress, trends and directions[END_REF], and the references therein).

Applying to System (1) the Smith predictor [START_REF] Smith | A controller to overcome dead time[END_REF], which is basic for controlling delay systems (see, e.g., [START_REF] Franklin | Feedback Control of Dynamic Systems[END_REF], [46]), yields:

1 ẏ(t + L) = ku(t) -d(t + L) (2)
where ŷ(t+L) and d(t+L) are respectively the estimates of y and d at time t + L. In order to obtain d(t + L), we need what may be viewed as a short-term forecast. Here it will be achieved via time series, as repeatedly in chain supply chain management (see, e.g., [START_REF] Axsäter | Inventory Control[END_REF], [START_REF] Aviv | A time-series framework for supply-chain inventory management[END_REF], [START_REF] Thomopoulos | Demand Forecasting for Inventory Control[END_REF]). We apply here, like [START_REF] Hamiche | A robust and easy approach for demand forecasting in supply chains[END_REF], a new setting which was proven to be useful in financial engineering [START_REF] Fliess | A mathematical proof of the existence of trends in financial time series[END_REF], [START_REF] Fliess | A-t-on vraiment besoin d'un modèle probabiliste en ingénierie financière[END_REF], and for motorway traffic [START_REF] Abouaïssa | On short-term traffic flow forecasting and its reliability[END_REF] and solar energy [START_REF] Fliess | Prediction bands for solar energy: New short-term time series forecasting techniques[END_REF] management. From a control-theoretic standpoint, this communication might be one of the very few publications that are combining Smith predictors, advanced forecasting techniques, and model-free control [START_REF] Fliess | Model-free control[END_REF] since Equation [START_REF] Abbou | On inventory control for perishable inventory systems subject to uncertainties on customer demands[END_REF], which is only a poor modeling, is mainly used for the purpose of computer simulations. Our result indicates how to mitigate the bullwhip effect at least in some cases. 2Several convincing computer simulations are displayed. Controlling perishable inventories3 is a not only a natural extension but also a key topic (see, e.g., [START_REF] Dolgui | Supply Chain Engineering[END_REF], [START_REF] Nahmias | Perishable Inventory Systems[END_REF], [START_REF] Sachs | Retail Analytics -Integrated Forecasting and Inventory Management for Perishable Products in Retailing[END_REF], and the references therein). In order to provide straightforward numerical simulations, replace according to [START_REF] Ignaciuk | Dead-time compensation in continuousreview perishable inventory systems with a remote supply[END_REF] System (1) by

ẏ(t) = -σy(t) + ku(t -L) -d(t), t ≥ L (3)
where σ 0 might be poorly known. The same tools as before are employed.

Our paper is organized as follows. Our main tools are briefly reviewed in Section II. Numerous computer simulations are presented in Section III. Some concluding remarks may be found in Section IV.

II. ABOUT OUR MAIN TOOLS A. Forecasting via time series

Our presentation differs quite a lot from the existing approaches, where other time series techniques are advocated (see, e.g., [START_REF] Axsäter | Inventory Control[END_REF], [START_REF] Bianchi | Recurrent Neural Networks for Short-Term Load Forecasting -An Overview and Comparative Analysis[END_REF], [START_REF] Brockwell | Time Series: Theory and Methods[END_REF], [START_REF] Mélard | Méthodes de prévision à court terme[END_REF], and the references therein).

1) Time series: Take the time interval [0, 1] ⊂ R and introduce as often in nonstandard analysis (see, e.g., [START_REF] Robinson | Non-standard Analysis[END_REF], [START_REF] Diener | Tutorial[END_REF], [START_REF] Lobry | Nonstandard analysis and representation of reality[END_REF], [START_REF] Fliess | A mathematical proof of the existence of trends in financial time series[END_REF]) the infinitesimal sampling

T = {0 = t 0 < t 1 < • • • < t ν = 1} where t i+1 -t i , 0 ≤ i < ν, is infinitesimal, i.e., "very small." A time series X(t) is a function X : T → R.
A time series X : T → R is said to be quickly fluctuating, or oscillating, if, and only if, the integral A X dm is infinitesimal, i.e., very small, for any appreciable interval, i.e., an interval which is neither very small nor very large.

According to a theorem due to Cartier and Perrin [START_REF] Cartier | Integration over finite sets[END_REF] the following additive decomposition holds for any time series X, which satisfies a weak integrability condition,

X(t) = E(X)(t) + X fluctuation (t) (4)
where

• the mean, or trend, E(X)(t) is "quite smooth,"

• X fluctuation (t) is quickly fluctuating. The decomposition (4) is unique up to an infinitesimal.

2) Forecasting: Let us start with the first degree polynomial time function p 1 (τ ) = a 0 + a 1 τ , τ ≥ 0, a 0 , a 1 ∈ R. Rewrite thanks to classic operational calculus with respect to the variable τ (see, e.g., [START_REF] Erdélyi | Operational Calculus and Generalized Functions[END_REF]) p 1 as P 1 = a0 s + a1 s 2 . Multiply both sides by s 2 :

s 2 P 1 = a 0 s + a 1 (5)
Take the derivative of both sides with respect to s, which corresponds in the time domain to the multiplication by -t:

s 2 dP 1 ds + 2sP 1 = a 0 (6)
The coefficients a 0 , a 1 are obtained via the triangular system of equations ( 5)-( 6). We get rid of the time derivatives, i.e., of sP 1 , s 2 P 1 , and s 2 dP1 ds , by multiplying both sides of Equations ( 1)-( 2) by s -n , n ≥ 2. The corresponding iterated time integrals are low pass filters which attenuate the corrupting noises. A quite short time window is sufficient for obtaining accurate values of a 0 , a 1 . Note that estimating a 0 yields the mean.

The extension to polynomial functions of higher degree is straightforward. For derivative estimates up to some finite order of a given smooth function f : [0, +∞) → R, take a suitable truncated Taylor expansion around a given time instant t 0 , and apply the previous computations. Resetting and utilizing sliding time windows permit to estimate derivatives of various orders at any sampled time instant. See [START_REF] Fliess | Non-linear estimation is easy[END_REF], [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF], [START_REF] Sira-Ramírez | Algebraic Identification and Estimation Methods in Feedback Control Systems[END_REF] for more details.

Set the following forecast X est (t + ∆T ), where ∆T > 0 is not too "large,"

X forecast (t + ∆T ) = E(X)(t) + dE(X)(t) dt e ∆T (7)
where E(X)(t) and dE(X)(t) dt e are estimated like a 0 and a 1 above. Let us stress that what we predict is the mean and not the quick fluctuations (see also [START_REF] Fliess | A mathematical proof of the existence of trends in financial time series[END_REF], [START_REF] Fliess | A-t-on vraiment besoin d'un modèle probabiliste en ingénierie financière[END_REF], [START_REF] Fliess | Prediction bands for solar energy: New short-term time series forecasting techniques[END_REF]).

B. Model-free control

The indisputable practical successes (see [START_REF] Abouaïssa | Energy saving for building heating via a simple and efficient model-free control design: First steps with computer simulations[END_REF], [START_REF] Bara | Toward a model-free feedback control synthesis for treating acute inflammation[END_REF], the numerous references therein, and [START_REF] Abouaïssa | On the control of robot manipulator: A model-free approach[END_REF], [START_REF] Afsi | Model-free control of a seeded batch crystallizer[END_REF], [START_REF] Barth | Fixed-wind UAV with transitioning flight capabilities: Model-Based or Model-Free Control approach? A preliminary study[END_REF], [START_REF] Bekcheva | Meilleure élasticité "nuagique" par commande sans-modèle[END_REF], [START_REF] Bekcheva | Cascaded model-free control for trajectory tracking of quadrotors[END_REF], [START_REF] Hong | A model-free control strategy for battery energy storage with an application to power accommodation[END_REF], [START_REF] Madadi | Comparison of different model-free control methods concerning real-ttme benchmark[END_REF], [START_REF] Menhour | An efficient model-free setting for longitudinal and lateral vehicle control: Validation through the interconnected Pro-SiVIC/RTMaps prototyping platform[END_REF], [START_REF] N'doye | Intelligent proportional-integral-derivative control-based modulating functions for laser beam pointing and stabilization[END_REF], [START_REF] Rampazzo | Modelling, simulation and real-time control of a laboratory tide generation system[END_REF], [START_REF] Ticherfatine | Fast ferry smoothing motion via intelligent PD controller[END_REF], [START_REF] Rocher | La production de nitrites lors de la dénitrification des eaux usées par biofiltration -stratégie de contrôle et de réduction des concentrations résiduelles[END_REF], [START_REF] Roman | Virtual reference feedback tuning of model-free control algorithms for servo systems[END_REF], [66], [START_REF] Telsang | Model-free control of building HVAC systems to accommodate solar photovoltaic energy[END_REF], [START_REF] Zhao | Double-objective finite control set model-free predictive control with DSVM for PMSM drives[END_REF], [START_REF] Zhou | Model-free deadbeat predictive current control of a surface-mounted permanent magnet synchronous motor drive system[END_REF]) of model-free control [START_REF] Fliess | Model-free control[END_REF] explains why it has already been summarized recently many times: [START_REF] Abouaïssa | On ramp metering: Towards a better understanding of ALINEA via model-free control[END_REF], [START_REF] Andréa-Novel | Model-free control of longitudinal and lateral dynamics for automated vehicles[END_REF], [START_REF] Bara | Toward a model-free feedback control synthesis for treating acute inflammation[END_REF], [START_REF] Fliess | Deux améliorations concurrentes des PID[END_REF], [START_REF] Menhour | An efficient model-free setting for longitudinal and lateral vehicle control: Validation through the interconnected Pro-SiVIC/RTMaps prototyping platform[END_REF]. It will therefore not be repeated here. This choice permits moreover to present computer experiments more thoroughly.

III. INVENTORY CONTROL SIMULATIONS

A. Classic situation 1) Model-based proportional feedback controller: The Smith predictor (2) stemming from Equation (1) yields the proportional controller [START_REF] Franklin | Feedback Control of Dynamic Systems[END_REF], or P controller,

u(t) = ẏ * (t + L) + d(t + L) -K p ê(t + L) k (8)
where y * (t) is a reference inventory trajectory and ê(t + L) = ŷ(t + L) -y * (t + L) is the predicted tracking error. This last quantity is obtained thanks to the forecast d(t + L) of the customer demand via the techniques sketched in Section II-A.2. Figure 1 displays excellent results with K p = 0.1. The data are borrowed from [START_REF] Schwartz | A process control approach to tactical inventory management in production-inventory systems[END_REF]. The customer demand is represented in Figure 1-(c): note that after the constant portion at the beginning a corrupting uniform white noise n, -1 ≤ n ≤ 1, has been added. The sampling time is 1 day. In Equation ( 1), L = 5, k = 0.95. With a poorly known model, i.e., k model = 0.95 × 0.9 = 0.855, a noticeable deterioration is shown in Figure 2 with the same assumptions, i.e., k = 0.95.

2) Model-free control and the corresponding intelligent controller: The obvious analogue of the ultra-local model in the model-free setting [START_REF] Fliess | Model-free control[END_REF] reads here:

ẏ(t) = αu(t -L) -F (t) (9)
where F (t) subsumes the poorly known internal system structure as well as the external perturbations. The coefficient α is chosen by the practitioner in such a way that the three terms in Equation ( 9) are of the same magnitude. It yields the corresponding intelligent proportional controller [START_REF] Fliess | Model-free control[END_REF], or iP, The forecast F (t + L) of F is again obtained via Section II-A.2. With the same k model = 0.95 × 0.9 = 0.855 and the same K p = 0.1, and with α = 1, the results depicted in Figure 3 improve a lot when compared to Figure 2.

u(t) = ẏ * (t + L) + F (t + L) -K p ê(t + L) α (10)
3) Bullwhip effect attenuation: Figure 3-(b) in particular shows the vanishing of the static error due to a poor knowledge of the model. According to the literature it should be viewed as an instance of bullwhip effect attenuation.

Assume now that the error due to the forecast is bounded by some quantity M > 0. Straightforward calculations show that the choice of a high gain (see, e.g., [46]) K p ≫ 0 minimizes the effect on y. Note however that the oscillations of the tracking error e are dangerously amplified on u by such a choice, which should be therefore avoided as much as possible. An "efficient" forecasting technique is therefore of utmost importance.

B. Perishable inventory

Consider System (3). Apply the iP [START_REF] Aviv | A time-series framework for supply-chain inventory management[END_REF]. Data are borrowed from [START_REF] Ignaciuk | Dead-time compensation in continuousreview perishable inventory systems with a remote supply[END_REF]. Choose L = 7, σ = 0.08, and k = 1. The sampling time is now 0.1 day.

Start with d(t) = 0 in Figures 4, 5, 6. The excellent results in Figure 4 lead to α = 1, K p = 0.1.

Inventory tracking is impressive. In order to check robustness, we selected respectively σ = 0.06, σ = 0.1 in Figures 5,6. The outcomes remain remarkable.

In Figure 7-(a) a demand with violent changes is introduced (see Fig. 7-(c)). The inventory tracking is still exceptional.

In Figure 8 a corrupting Gaussian white noise, with a standard deviation equal to 10, is taken into account (see Fig. 8-(c)). In order to avoid negative values for u(t), d(t) oscillates around the value 30. The results are still convincing.

Remark 1: Similar results to those in Section III-A.3 on the bullwhip effect could be extended here without any difficulty.

IV. CONCLUSION

A true application to concrete situations is necessary of course in order to confirm the above approach. As in financial engineering risk analysis would also be a stimulating topic in inventory control [START_REF] Schlegel | Supply Chain Risk Management -An Emerging Discipline[END_REF]: the tools developed in [START_REF] Fliess | Prediction bands for solar energy: New short-term time series forecasting techniques[END_REF], [START_REF] Fliess | A-t-on vraiment besoin d'un modèle probabiliste en ingénierie financière[END_REF] might lead to a solution. Related means [START_REF] Fliess | Prediction bands for solar energy: New short-term time series forecasting techniques[END_REF] should also lead to appropriate metrics (see, e.g., [START_REF] Trapero | A novel time-varying bullwhip effect metric: An application to promotional sales[END_REF] and the references therein).

From a control-theoretic viewpoint, a completely new approach combining Smith predictors, short-term forecasting and model-free control has been presented in order to provide an adequate regulation of some delay systems. Let us remind [START_REF] Fliess | Model-free control[END_REF] the reader about the severe difficulties that are due to delays in the model-free setting. Compare with other techniques in the literature: [START_REF] Ahmed | Model-free control using time delay estimation and fractional-order nonsingular fast terminal sliding mode for uncertain lower-limb exoskeleton[END_REF], [START_REF] Join | Vers une commande sans modèle pour aménagements hydroélectriques en cascade. 6 e Conf. Int. Francoph. Automat[END_REF], [START_REF] Doublet | Model-free control for unknown delayed systems[END_REF], [START_REF] Doublet | Stability analysis for unknown delayed systems controlled by model-free control[END_REF], [START_REF] Han | Model-free based adaptive nonsingular fast terminal sliding mode control with time-delay estimation for a 12 DOF multi-functional lower limb exoskeleton[END_REF], [START_REF] Huba | Tuning of an iP controller for a time-delayed plant by a performance portrait method[END_REF], [START_REF] Ghazally | Vibration control of an active vehicle suspension systems using optimized model-free fuzzy logic controller based on time delay estimation[END_REF], [START_REF] Thabet | Design of adaptive PID controllers based on adaptive Smith predictor for ultra-local model control[END_REF], [START_REF] Yaseen | Synthesis of model-free control for system with time-varying communication delay[END_REF], [START_REF] Yaseen | A model-free approach to networked control system with time-varying communication delay[END_REF], [START_REF] Zhang | Model-free based neural network control with time-delay estimation for lower extremity exoskeleton[END_REF].

Fig. 2 :

 2 Fig. 1: Scenario 1

Fig. 3 :Fig. 4 :Fig. 5 :Fig. 6 :Fig. 7 :Fig. 8 :

 345678 Fig. 3: Scenario 3

Smith predictors were already employed a few times for supply chain management. See, e.g.,[START_REF] Chołdowicz | Comparison of a perpetual and PD inventory control system with Smith predictor and different shipping delays using bicritical optimization and SPEA2[END_REF],[START_REF] Ignaciuk | Smith predictor based control of continuous-review perishable inventory systems with multiple supply alternatives[END_REF],[START_REF] Ignaciuk | Dead-time compensation in continuousreview perishable inventory systems with a remote supply[END_REF],[START_REF] Escudero-Gomez | Inventory trajectory tracking via an demand-rate-estimation using Smith's predictor[END_REF].

The literature does not seem to contain any clear-cut definition of this effect.

Fresh food for instance.


Franc ¸ois Haibo Zhang

Diana M Gruson

Christophe Flórez Rodriguez

Saudemont

Overvoltage Limitation Method of an Offshore Wind Farm with DC Series-Parallel Collection Grid

Keywords: offshore wind farm, HVDC transmission, modular multilevel converter (MMC), overvoltage limitation

This paper describes the characteristics of a series parallel wind farm (SPWF) topology and investigates the control strategy to ensure its safe operation. The SPWF was found to have advantages over other pure dc wind farm architectures in terms of lower construction cost and lower power losses in the collection system. However, unbalance power productions among the wind turbines cause the variations of their output voltages, which may endanger the safe operation of the entire wind farm. This paper proposes a global control strategy that prevents wind turbines from operating above their overvoltage capabilities. With an active participation of the onshore converter, the proposed strategy allows maximum power point tracking (MPPT) of the wind turbines. The practical limitations of this strategy are discussed and improvements are given. The feasibility of the proposed method is validated in a simulation of 300 MW wind farm developed in EMTP-RV.

I. INTRODUCTION

T HE exploitation of far-offshore wind energy introduces new challenges in offshore power transmission and ac grid integration. It has been found unsuitable to use HVAC for long distance submarine transmission, because the reactive power generated by ac cables can considerably undermine its transmission capability [START_REF] Negra | Loss evaluation of HVAC and HVDC transmission solutions for large offshore wind farms[END_REF]. Consequently, high voltage dc (HVDC) technology is dominating in large long-distance offshore projects.

Conventional point-to-point HVDC transmission increases the power conversion stages compared to HVAC, which indicates higher power conversion losses and higher installation and maintenance costs [START_REF] Bresesti | Hvdc connection of offshore wind farms to the transmission system[END_REF]. Correspondingly, potential efficiency improvement and cost reduction reside in eliminating large power transformers, power converters and platforms. This leads to recent investigations on pure dc offshore wind farm with dc collection and transmission systems [START_REF] Meyer | Control and design of dc grids for offshore wind farms[END_REF], [START_REF] Lundberg | Wind farm configuration and energy efficiency studiesseries DC versus AC layouts[END_REF].

Most pure dc configurations simply replace the power transformers in the wind turbine and in the offshore substation with dc/dc converters [START_REF] Max | System efficiency of a dc/dc converter-based wind farm[END_REF]. These configurations do not reduce the power conversion stages. Instead, the focus is placed on improving the power conversion efficiency by new dc/dc converter designs [START_REF] Holtsmark | An all-dc offshore wind farm with series-connected turbines: An alternative to the classical parallel ac model?[END_REF]- [START_REF] Chen | Analysis and comparison of medium voltage high power dc/dc converters for offshore wind energy systems[END_REF].

Through recent studies, a promising solution to reduce the power conversion stages of the pure dc configurations is to series connect wind turbines with dc output in a cluster. The cluster establishes the HVDC transmission voltage directly at the collection grid side, so that the centralized offshore converter and its platform are not necessary. The wind turbines in the cluster must share a same cluster current, which gives rise to a group of series topologies using thyristor based HVDC [START_REF] Nishikata | A new interconnecting method for wind turbine/generators in a wind farm and basic performances of the integrated system[END_REF]. A HVDC link using a current source inverter (CSI) at onshore and diode rectifiers at offshore is studied in [START_REF] Kawaguchi | Offshorewind-farm configuration using diode rectifier with mers in current link topology[END_REF]. A similar thyristor based configuration using simple buck choppers at the wind turbine outputs is proposed in [START_REF] Veilleux | Interconnection of direct-drive wind turbines using a series-connected DC grid[END_REF]. The main disadvantages of these configurations are the limitations of line commutated thyristors [START_REF] Flourentzou | Vsc-based hvdc power transmission systems: An overview[END_REF], and the usage of transformer-less converters in the wind turbines, so that there is no galvanic isolation between the wind turbine and the HVDC transmission system.

Another potential method to series connect the wind turbines is to use voltage source converter (VSC) based HVDC configuration, which is firstly proposed in [START_REF] Lundberg | Wind farm configuration and energy efficiency studiesseries DC versus AC layouts[END_REF] called series parallel wind farm (SPWF) topology, and this is the method considered in this paper. The studied power generation and transmission system is shown in Fig. 1. The HVDC voltage is controlled by an onshore multilevel modular converter (MMC). In the collection grid, full bridge dc/dc converters are used to step up the low wind turbine generator side dc voltage to medium level. In this paper, it is assumed that the intermediate transformer in the dc/dc converter also serves to provide galvanic isolation for the wind turbine. However, the galvanic isolation level of the first and last wind turbines to the common mode voltage, equivalent to half of the HVDC transmission voltage, still remains a practical challenge to the academic and industrial communities. This paper focuses on the operation issues of this MMC-HVDC based SPWF, aiming to increase the power availability of this topology.

The comparison of the SPWF with other dc configurations regarding power losses, costs of components, and power availability/reliability are studied in [START_REF] Holtsmark | An all-dc offshore wind farm with series-connected turbines: An alternative to the classical parallel ac model?[END_REF], [START_REF] Bahirat | Comparison of wind farm topologies for offshore applications[END_REF]. The clear advantages of the SPWF are lower collection grid losses, and lower construction and maintenance costs, thanks to the elimination of the offshore converter and platform. However, the MMCcontrolled HVDC voltage imposes overvoltage on certain series connected wind turbines operating in uneven wind conditions. If the overvoltage surpasses a permitted limit, it may lead to cascade failure of the entire cluster. A solution is suggested in [START_REF] Lundberg | Wind farm configuration and energy efficiency studiesseries DC versus AC layouts[END_REF] to utilize an output voltage limitation controller which actively reduces the power production of overvoltage units. However, this leads to undesired power losses. Another solution proposed in [START_REF] Chuangpishit | Topology design for collector systems of offshore wind farms with pure DC power systems[END_REF] can mitigate the overvoltage level by topology modification. The modified collection grid allows currents to flow into adjacent clusters. But the substantial investment on extra switch gears and overrated cables is not considered in the author's cost analysis. This paper, therefore, proposes a novel overvoltage limitation strategy which requires communication between wind farm and onshore converter, and active participation of the onshore converter in limiting the wind turbines output voltages. Section II of this paper describes the studied system. The proposed limitation strategy, practical limitations and improvements are presented in section III. Section IV deals with the corresponding control systems. In Section V, a simulation on a 300 MW wind farm has been developed in EMTP-RV [START_REF] Mahseredjian | On a new approach for the simulation of transients in power systems[END_REF] to investigate the feasibility of proposed control strategy. The conclusion is made in Section VI.

i la . . . R arm L arm R arm L arm . . . i ua i HVDC SM Nla SM 2la SM 1la SM Nua SM 2ua SM 1ua i lb . . . R arm L arm R arm L arm . . . i ub SM Nlb SM 2lb SM 1lb SM Nub SM 2ub SM 1ub i lc . . . R arm L arm R arm L arm . . . i uc SM Nlc SM 2lc SM 1lc SM Nuc SM 2uc SM 1uc i gc i gb N R g L g i ga + - v ga C SM + - u Ci + - u SM PMSG DC AC DC DC L SR WT m,1 + - u out,m,1 PMSG DC AC DC DC L SR WT 3,1 + - u out,3,1 PMSG DC AC DC DC L SR WT 2,1 + - u out,2,1 PMSG DC AC DC DC L SR WT 1,1 + - u out,1,1 i CLn i CL2 i CL1 Cluster 2 . . . Cluster n + - u HVDC . . .

II. SYSTEM DESCRIPTION

Fig. 1 shows the entire studied system in which a m × n SPWF is integrated to the onshore grid by MMC based HVDC technology. The collection grid of the SPWF is composed of n parallel connected clusters. In each cluster, m wind turbines are connected in series to establish the HVDC voltage. The MMC consists of six arms of cascaded submodules (SMs). A SM with half bridge topology is shown in the figure.

The whole system features the modular multilevel concept at both offshore and onshore sides. But there is substantial difference between the SPWF and MMC from the energy distribution point of view. Energy balancing among the SM capacitors is one of the main problems of MMC. The balancing capacitor algorithm (BCA) is implemented, so that the SM capacitor voltages are regulated equally [START_REF] Tu | Reduced switching-frequency modulation and circulating current suppression for modular multilevel converters[END_REF]. However, the power of each wind turbine in the SPWF depends on the wind speed, which is not controllable. Unequal power production among the units leads to uneven allocation of output voltage across the dc/dc converters. Therefore, the output dc/dc converters are required of variable output voltage capability, and thus, variable input-output transformation ratio. The significance of this characteristics will be further addressed.

A. SPWF main components

The wind turbines employ permanent magnetic synchronous generators (PMSGs) and fully rated power rectifiers. The dc/dc converter steps up the low intermediate dc link voltage to medium voltage. The most adopted dc/dc converters with insulation transformers in the literature are: full bridge converter (FBC) [START_REF] Sabate | Design considerations for high-voltage high-power full-bridge zero-voltageswitched pwm converter[END_REF], single active bridge (SAB) [START_REF] De Doncker | A threephase soft-switched high-power-density dc/dc converter for high-power applications[END_REF], dual active bridge converter (DAB) [START_REF] Kheraluwala | Performance characterization of a high-power dual active bridge dcto-dc converter[END_REF] and series resonant converter (SRC) [START_REF] Steigerwald | A comparison of half-bridge resonant converter topologies[END_REF]. In this paper, the FBC is adopted for two reasons: 1) the buck type converter ensures a continuous flow of the cluster current; 2) the ability to vary transformation ratio. The FBC is shown in Fig. 2a and corresponding duty cycle control block is shown in Fig. 2b. The outer loop of the control aims to maintain the intermediate dc link voltage u in to be constant, while the output voltage of the FBC is not controlled.

To prevent the failure of a single unit from triggering cascade failure in the cluster, a bypass protection system is needed to bypass the fault wind turbine. Two bypass solutions are therefore proposed. The first solution is shown in Fig. 3a. The dc breakers at the output of the dc/dc converter directly cut off the entire wind turbine from the cluster. Simultaneously, a connector switches on to provide a current path. The second solution is to use the diode bridge of the FBC as the cluster current path, as shown in Fig. 3b. When the internal fault is detected, the IGBTs are kept at off-state. The output capacitor voltage gradually decreases to zero.

The smoothing reactor in the cluster is designed to a value big enough to limit the gradient of current boost ∆i HVDC /∆t during critical HVDC faults. The smoothing reactor is split evenly to as many parts as the number of wind turbines, and thus each part L SR is small enough to be placed inside the wind turbine tower, as shown in Fig. 1.

i in + - u in i pc i p + - u p + - u s i out L out C out i CL + - u sc + - u out 1 : σ (a) Converter circuit K puin + K iuin s K piout + K iiout s u * in + u in - - i in + i * pc × σ•D ÷ i * out + i out - + u out + u * sc × σ•u in ÷ D (b) Dual-loop control blocks

B. Steady state operation principles

In the collection system, the HVDC link voltage is regulated by the onshore MMC, while the wind turbine output voltage is not individually controlled. Assuming that the power electronics devices in wind turbines are ideal and have infinite voltage and current capabilities, the operation of wind turbines follows the Kirchhoff Circuit Laws. Output voltage of the x th wind turbine located in cluster y, WT x,y , is calculated as:

u out,x,y = u HVDC p out,x,y m x=1 p out,x,y = u HVDC p out,x,y p CL,y (1)
i CL,y = p CL,y u HVDC = p out,x,y u out,x,y (2)
where the wind turbine power production and output voltage are represented by p out,x,y and u out,x,y . p CL,y and i CL,y are the power production and current of cluster y, respectively. Equation ( 1) gives the steady state output voltage of a wind turbine in the SPWF. It implies that the output voltage depends not only on its own power, but also on the power of its cluster. It has been mentioned that the output dc/dc converter regulates its input dc link voltage instead of its output voltage. The output voltage variation is illustrated in Fig. 4. When the wind turbine power generation is higher or lower than the average power in the cluster, the power unbalance is reflected by the unmatched currents at both sides of the wind turbine output capacitor. The stiff current i out charges or discharges the output capacitor, which consequently results in the increase or decrease of the output voltage u out , as shown in Fig. 4.

C. Definition of overvoltage

The HVDC voltage is regulated to m times the wind turbine nominal voltage U out,nom by the onshore MMC. Due to the wind speed variation, certain wind turbines output voltage may rise above its nominal value, and these particular units are in overvoltage. Bypassing fault wind turbines also pushes the remaining fault-free units to operate in overvoltage. The maximum voltage level which a wind turbine converter can support is its overvoltage capability U limit . The ratio between the overvoltage capability and the nominal output voltage is denoted as α.

A higher overvoltage capability permits a larger power disparity in the cluster. Nevertheless, physically expanding the size of converter is neither economical nor practical. For these reasons, overvoltage limitation control should be developed to prevent wind turbine from operating above their permitted overvoltage capability.

III. OVERVOLTAGE LIMITATION CONTROL

As mentioned in the introduction, a local overvoltage limitation method is proposed in [START_REF] Lundberg | Wind farm configuration and energy efficiency studiesseries DC versus AC layouts[END_REF]. When the output voltage reaches the overvoltage capability, an output voltage limitation controller is triggered to prevent the output voltage from increasing, by the means of reducing the wind turbine power production. The control diagram is shown in Fig. 5. The output of the voltage limitation controller is taken as a compensation torque T em,com to the inner current controller of the rectifier. The value of this torque will always be negative, and thus the torque reference for the inner current loop reduces. The compensation torque disables the MPPT of the wind turbine, leading to wind power curtailment.

In this paper, a new solution named as global overvoltage limitation strategy is proposed, which ensures both the safe

A. Basics of the global overvoltage limitation strategy

The proposed global control strategy divides the wind farm operation into two working modes: normal mode and overvoltage limitation mode. In normal mode, the wind turbines output voltages are allowed to vary under their overvoltage capabilities U limit . The reference HVDC voltage u * HVDC of the MMC is constant, equivalent to the sum of wind turbines nominal output voltages U out,nom :

u * HVDC = mU out,nom (3)
When the output voltage of any wind turbine reaches its overvoltage capability level, the onshore converter enters into voltage limitation mode. The HVDC voltage is regulated to the level which can ensure the output voltage of the wind turbine with the maximum power production in the wind farm to stay under the the overvoltage limit:

u * HVDC = min U out,nom p CL,y max(p out,x,y ) . ( 4
)
The control scheme of the global control and corresponding controller references of the MMC are depicted in Fig. 6.

The bypassed fault unit can be regarded as a wind turbine with zero power production and therefore its output voltage equals to zero. The corresponding cluster is reduced to a series connection with (mk) units, k is number of bypassed units. In this case, all the remaining units in the cluster are operated at overvoltage, if the HVDC voltage remains at mU out,nom . Even all the rest units in the cluster produce nominal power, their output voltages are clamped to m m-k U out,nom . With this clamped voltage as base voltage, there is very little margin

u base HVDC = m -k m U out,nom . (5)
The overvoltage limitation strategy obeys the above working modes based on this reduced HVDC voltage.

The advantage of this strategy is that wind turbines do not need to reduce their own power production. Accurate information of all wind turbines operating condition is essential, which depends on data communication via the fibre optic cable embedded in the HVDC cables.

B. Constraints of basic global strategy

It is important to note some practical limitations of the proposed global control strategy:

1

) Communication delay and interruption

The communication delay is of the order of several dozens of milliseconds. In large power systems such as in the Bonneville Power Administration system, the latency of fibre communication is reported as approximately 38 ms for one way [START_REF] Wu | Evaluation of time delay effects to wide-area power system stabilizer design[END_REF]. The latency is caused by the fact that the electrical devices do not have inherent communication capabilities, and therefore they have to rely on embedded computer systems to serve as communication interfaces. The most significant delays during the process are the data acquisition delay, the packet processing and transmission delay [START_REF] Wang | A survey on the communication architectures in smart grid[END_REF]. In the SPWF, the MMC response time and capacitive HVDC cables also lead to the event responding delay.

2) Minimum allowable HVDC voltage level

The second constraint is related to the practical limitations of the MMC. In order to prevent the diode bridge from forward bias and consequently lose the full control of the MMC, the HVDC voltage u HVDC should always be greater than the peak value of ac side phase-to-phase voltage:

u * HVDC > √ 6V g (6)
where V g is the RMS value of the grid phase-to-neutral voltage.

3) HVDC cables and onshore MMC losses

The proposed global control reduces the HVDC voltage, which definitely increases the HVDC cable current and cable power losses. Increased HVDC current also has an influence on the conduction and switching losses of the semiconductors in the MMC. The MMC has a much lower switching frequency compared to a conventional 2-level or 3-level VSC. It is found in [START_REF] Freytes | Losses estimation method by simulation for the modular multilevel converter[END_REF] that the biggest share of power losses presented in the MMC is the conduction losses of the semiconductors. Usually, one suggests that it is necessary to limit the RMS value of the currents flowing through the converters to limit their losses. However, reference [START_REF] Gruson | Impact of control algorithm solutions on modular multilevel converters electrical waveforms and losses[END_REF] shows that the conduction losses of the MMC are not only related to the currents' RMS values, but also the waveforms of these currents and therefore the MMC control algorithms. Because the total losses in the MMC are very low (around 1% of the nominal power), the estimation of MMC losses is not studied in this paper.

C. Optimal HVDC voltage

Due to the above issues of the global control, it is suggested that the local overvoltage limitation controllers be installed in the wind turbines. When the global control reaches the mentioned limitations, the local control should be activated to serve as complementary limitation controllers. To summarize, the two kinds of overvoltage limitation controls play different roles in the SPWF:

• The local overvoltage limitation control activates during the communication delay and communication interrupts. Furthermore, it serves as a back-up control for when the HVDC has been reduced to its minimum allowable level.

• The global control serves to limit the wind turbines output voltages without power curtailment by reducing the HVDC voltage. It is selected as the primary option for voltage limitation. The process of controlling the series wind farm by these two methods is illustrated in Figure 7.

The proposed global control increases the power production of the wind farm, comparing to the local control. However, the HVDC cable losses are also increased. This leads to an investigation of the optimal HVDC voltage level, which enables the onshore side to receive maximum power from the offshore side:

• On one hand, if only local control strategy is applied, this situation causes the maximum wind turbine power curtailment and the minimum cable losses. • On the other hand, if only the global control is applied, this situation causes the maximum cable losses and the minimum wind turbine power curtailment. As a result, the optimal voltage ûopt HVDC has a compromised level, which must drop in the range:

ûgctl HVDC ≤ ûopt HVDC ≤ U HVDC,nom (7)
where the hat ˆdenotes the variables at a reduced HVDC voltage level. ûgctl HVDC is the reduced HVDC voltage if only the global strategy is applied, U HVDC,nom is the nominal HVDC voltage.

In order to find this optimal voltage level, we decrease the HVDC voltage gradually, and then calculate the instant wind farm power production pWF at each HVDC voltage level (note that the wind farm power production is partially reduced due to the accompanied local limitation controllers), as well as the cables losses pcable . At each voltage level, the subtraction of the wind farm power production and the cables power losses is the amount of power received at onshore:

preceived = pWF -pcable (8)

D. Case study

To illustrate the above method to find the optimal HVDC voltage, a 20 × 1 SPWF is used for example. The nominal wind turbine output power P out,nom and voltage U out,nom are 5 MW and 32 kV, overvoltage ratio α = 0.1. For a practical comparison of the cables power losses and the power curtailment due to local control, the commissioning Dolwin1 project is taken as reference. The power transmission length is 165 km. HVDC transmission voltage is ±320 kV. The cables with copper conductor should have a section area of 85 mm 2 according to [START_REF] Abb | HVDC light -cables submarine and land power cables[END_REF].

Assume that:

U out,nom = P out,nom = I CL,nom = 1 pu; U limit = 1.1 pu
The SPWF power instant production and output voltages without control strategy applied are:

p out = WTx,1           1-4
1.00 Power received at onshore end Fig. 8: Variation of the cable losses, SPWF generated power, and transmitted power to onshore. received power at every HVDC voltage level, the variation of the received power is plotted in Figure 8.

At nominal HVDC voltage 640 kV, both the cable losses and the generated power from the SPWF are minimal. The local control needs to reduce the power production of WT 1-8,1 from 1 pu to 0.9036 pu, equivalent to 3.85 MW power curtailment. The power received at onshore is 80.8 MW. As expected, the cable losses increase as the HVDC voltage decreases. The received power peaks at 604.8 kV to 84.3 MW, 3.5 MW greater than using only local control strategy. Note that this voltage equals to the ûgctl HVDC where the local control is fully deactivated. This indicates that in this case, in the whole range of [û gctl HVDC , U HVDC,nom ], the reduction of voltage offers greater power generated by the SPWF than the cable losses. Further reducing the HVDC voltage below 604.8 kV, the cable losses increase while all wind turbines operate at MPPT, and therefore the power received at onshore side declines.

IV. SPWF-MMC MODELS AND CONTROL SYSTEMS

Simulating the detailed EMTP-RV models of the wind turbine and onshore MMC is time consuming. Without losing generality, simplified models are adopted in the simulation: 1) a simplified wind turbine model with basic wind power conversion, dc/dc converter is fundamental to study the operating principles and characteristics of a SPWF with a number of units; 2) a simplified MMC model with dynamics and stored energy in each arm is necessary to realize the global control.

A. Wind turbine generation system

This paper assumes that in normal mode, the PMSG wind turbine is controlled ideally to operate at MPPT. Based on this assumption, a simplified series wind turbine model can be created as shown in Fig. 9. The dc link voltage is controlled to be constant by the dc/dc converter, and thus the variation of the

i in i pc + - u in - + L out i out C out + - u out i in i pc + - u in - + L out i out C out + - u out Fig. 2b u in K p lctl + K i lctl s p em × u in ÷ p w + + U limit - u out +
Fig. 9: Simplified wind turbine model local overvoltage limitation controller. electromagnetic power p em is reflected as current variation of the controlled current source i in . The dc/dc converter control follows the dual loop control described in Fig. 2b. The local overvoltage limitation control is represented by an negative input power compensation. In normal operation, the output voltage is lower than the overvoltage limitation, and thus the input to the output voltage limitation PI controller is always 0. As the output voltage reaches the limitation, this controller is activated. A negative compensation is added to the maximum extractable power p w , and therefore p em decreases.

i la - + u la R arm L arm R arm L arm - + u ua i ua i HVDC i lb - + R arm L arm R arm L arm - + i ub i lc - + R arm L arm R arm L arm - + i uc i Ctotlc + - u Ctotlc C tot + - u Ctotuc i Ctotuc C tot m la m lb m lc m ua m ub m uc i gc i gb R g L g i ga + - v ga

B. Multilevel modular converter

In the literature, there exist several types of MMC models. This paper adopts the MMC arm average model, which assumes that all SM capacitor voltages are balanced. Although in this model the information of the switching states is lost, the dynamics and stored energy in each arm are retained [START_REF] Antonopoulos | On dynamics and voltage control of the modular multilevel converter[END_REF]. Correspondingly, the low level control of MMC SMs is not considered. Only the high level control regarding to the MMC currents, voltages and stored energy is given.

The N+1 level detailed MMC shown in Fig. 1 with n the number of active cells in the arm, new variables are deduced:

i * g d+ + v g d i * gq + + Q * g i g d i gq - - + + v g d + v gq + u * v d u * vq θ Arm energy controller ÷ i dif f controller u 2 * ctotj + u 2 ctotj - P * DCj i * dif fj + i dif fj
u ulj = m ulj u Ctot ulj i Ctot ulj = m ulj i ulj (9)
where the subscripts u and l represent the upper arm and the lower arm, and j referring to three phases a, b, c. By performing the following variables changes,

u dif fj = u uj + u lj 2 u vj = u lj -u uj 2 i dif fj = i uj + i lj 2 L eq = L g + 1 2 L arm R eq = R g + 1 2 R arm
the system can be decoupled to a dc part and an ac part (in dq reference frame) as:

u HVDC 2 -u dif fj = L arm di dif fj dt + R arm i dif fj (10)
u v d -v g d = L eq di g d dt + R eq i g d -ω g L eq i gq (11)
u vqv gq = L eq di gq dt + R eq i gq + ω g L eq i g d

Neglecting the energy stored in the inductances, power exchange between the dc side and the ac side of the MMC in one phase induces the instantaneous energy change in one MMC phase (upper arm and lower arm):

p dcj -p acj = dW Σ j dt = 1 2 C tot du 2 Ctotuj + du 2 Ctot lj dt . (13)
Considering the HVDC cables as small capacitance:

1 2 C dc du 2 HVDC dt = p WF -p dc . (14)
The entire MMC control system is shown in Fig. 11. Thanks to the decoupling of the MMC system, the control of MMC can be separated into inner current loop and outer power (voltage) loop. This part of control system is similar to conventional 2-level VSC, and therefore will not be addressed in detail. The outer loop controllers are tuned 10 times slower than the inner loop. An unique feature of MMC is that the energy is distributed in the MMC arms. As a result, an extra arm energy control loop is added as shown in Fig. 11.

In the traditional two-level or three-level VSC, the dc bus voltage is directly established by the energy stored in dc link capacitor. In contrast, the energy stored in the arms of MMC is not directly linked to the dc bus voltage. By setting the energy reference of u 2

Ctot ulj to the square of the dc bus voltage u 2 HVDC , the MMC stored energy can be related to the dc bus as studied in [START_REF] Samimi | MMC stored energy participation to the DC bus voltage control in an HVDC link[END_REF]. Equation ( 14) becomes:

dW Σ dt = 1 2 C tot 6 du 2 HVDC dt = p dc -p ac (15)
V. SIMULATION RESULTS

The simulation of a 20 × 3 dc SPWF with the proposed control strategy and corresponding MMC control is fully developed and implemented in EMTP-RV. For the wind turbines, the nominal power and output voltage are 5 MW and 32 kV. The nominal HVDC voltage for the MMC is 640 kV. The overvoltage ratio α = 0.1. The HVDC cables used in the simulation are 100 km WideBand models available in the EMTP-RV library. More parameters of SPWF, MMC and control systems are summarized in Tab. I and Tab. II.

A. Wind speed variation

The first scenario simulates the functionalities of the overvoltage limitation control strategies when wind speed varies. The maximum extractable wind powers p w are plotted in Figure 12a. Wind turbines in the cluster 2 and 3 have identical wind speed variations as the units in the cluster 1. All units operate under their overvoltage capabilities until 4.6 s, as the power production unbalance increases, WT 1-4,y reach their overvoltage capabilities.

In order to show the different behaviors of the local limitation control and the global limitation control, the communication delay is intentionally extended. Only the local limitation control is activated before 8 s and the global limitation is applied during 8 s to 10 s.

During 4.6 s to 8 s, the local limitation controllers in WT 1-4,y reduce their power production to 4.65 MW as shown in Figure 12b. The local limitation controllers succeed in regulating the output voltages of WT 1-4,y at 1.1 pu. The power received at onshore converter is -252.5 MW. The negative sign indicates that the power is transmitted from the HVDC system to the ac grid. From 8 s, the global limitation control is activated. The HVDC voltage is reduced to the optimal HVDC voltage ûgctl HVDC 604.8 kV. The power productions of WT 1-4,y are restored to nominal power production because the local control is released owing to the reduction of HVDC voltage. As shown in Figure 12d, the received power at onshore steps up to -256.7 MW. Comparing these two overvoltage limitation strategies, 4.2 MW of power is recovered thanks to the proposed method.

The HVDC voltage and current are plotted in Fig. 12e and12f. As explained in Subsection IV-B that, the energy reference of u 2

Ctot ulj is set to the square of the dc bus voltage u 2 HVDC , and therefore the MMC stored energy can be associated with the HVDC link. The variation of the MMC arm stored energy during the HVDC voltage reduction is plotted in Fig. 13. It can be noticed that the above description is validated.

B. One fault unit bypassed in the cluster 2

In the second scenario all the wind turbines in the wind farm produce nominal power. One fault unit is bypassed in cluster 2 at 5 s. The communication delay is not intentionally increased.

Fig. 14a only shows the power variation of the entire cluster 1 (same as cluster 3) and cluster 2. At 5 s, cluster 2 suffers a sudden drop of power since one unit is bypassed. At almost the same time, the MMC regulates the HVDC voltage to 19 pu in order to remove the overvoltage imposed on all remaining wind turbines in the cluster 2. This results in a decrease of output voltage in the units located in cluster 1. Because the wind turbine output capacitors and the smoothing reactors construct a resonance circuit, the power unbalance between parallel connected clusters causes voltage oscillations across each unit. Fig. 14c presents the output voltage of WT 1,1 and WT 1,2 . However, these oscillations do not occur on the transmission system. Fig. 14b and 14d plot the voltage and current of the HVDC transmission system. The energy only exchanges between the output capacitors and the smoothing reactors in the collection grid.

VI. CONCLUSION

A dc collection grid with series connection of wind turbines allows a reduction of the power conversion stages, and consequently construction cost and power losses in the collection grid. A dc/dc converter with variable transformation ratio and a fault bypass system are essential for a wind turbine to operate in the series topology which is connected to a VSC based HVDC system. Unbalanced power production among wind turbines will cause variation of their output voltages, and may even render certain wind turbines unacceptable overvoltage. This paper compares two output overvoltage limitation methods. A local control method is fast but leads to wind power curtailment. On the contrary, the proposed global control method can ensure MPPT control for all wind turbines, but has some practical limitations regarding communication delay and increased cable losses. It is suggested in this paper to combine both control methods to ensure safe operation of all wind turbines. And it is found that the power gained by using the global control overwhelms the increased power losses in the cables. A series parallel offshore wind farm of 300 MW with 60 wind turbines is simulated in EMTP-RV, which validates the proposed control method by showing accordance results in the collection system and in the onshore MMC with the analysis.

APPENDIX

  H. Zhang, F. Gruson, Diana M. Flórez R. and C. Saudemont are with Univ. Lille, Arts et Metiers ParisTech, Centrale Lille, HEI, EA 2697 -L2EP -Laboratoire d'Electrotechnique et d'Electronique de Puissance, F-59000 Lille, France.

Fig. 1 :

 1 Fig. 1: MMC based dc SP offshore wind energy generation system.

Fig. 2 :Fig. 3 :

 23 Fig. 2: Full bridge DC/DC converter and control blocks.

iFig. 4 :

 4 Fig.4: Variation of output voltage due to power imbalance between the wind turbine and average level of the cluster.

Fig. 5 :

 5 Fig. 5: Local overvoltage limitation control diagram.

Fig. 6 :

 6 Fig. 6: Global control scheme of the SPWF-MMC HVDC.

Fig. 7 :

 7 Fig. 7: Arrangement of global and local controls of the SPWF.

  to (1), WT 1-8,1 produce nominal power while their output voltage are greater than their capabilities 1.1 pu. By decreasing the HVDC voltage gradually and calculating the

Fig. 10 :

 10 Fig. 10: MMC arm average model.

  can be reduced to the average model shown in Fig. 10. The equivalent MMC arm capacitance, referring to C tot , equals to the sum of all SM capacitors C/N . By defining the duty cycle m = n/N

vFig. 11 :

 11 Fig. 11: High level control structure of MMC for HVDC voltage control.

Fig. 12 :

 12 Fig. 12: Simulation results of the SPWF with local and global overvoltage limitation methods.

2 CtotcFig. 13 :

 213 Fig. 13: MMC arm stored energy variation.

  Output voltage of wind turbine WT1,1 and WT1,2

Fig. 14 :

 14 Fig. 14: Simulation results of the SPWF when a fault unit is bypassed.

TABLE I :

 I System parameters

	Parameter	Notation	Value
	DC link capacitor	C in	10 mF
	Output capacitor	Cout	293 µF
	Output inductor	Lout	100 mH
	Smoothing reactor per WT	L SR	5 mH
	MMC arm resistor	Rarm	1.02 Ω
	MMC arm inductor	Larm	50 mH
	Grid resistor	Rg	50 mΩ
	Grid inductor	Lg	60 mH
	Arm equivalent capacitor	Ctot	25 µF

TABLE II :

 II Controllers parameters

	Parameter	Kp	K i
	dc/dc current loop	32	9000
	dc/dc voltage loop	0.42	9
	Local voltage limitation loop PLL Grid current loop Differential current loop	7300 100 20	1.8 × 10 5 3.6 × 10 5 2 × 10 4 4500
	HVDC voltage loop	0.004	0.08
	Stored energy loop	0.001	0.023


Dr Juan Manuel Asensio
email: asensior@insa-toulouse.fr

Julien Marbaix

Nicolas Mille

Lise-Marie Lacroix

Katerina Soulantica

Pier-Francesco Fazzini

Julian Carrey

To Heat or not to Heat: a Study of the Performances of Iron Carbide Nanoparticles in Magnetic Heating

Heating magnetic nanoparticles with high frequency magnetic fields is a topic of interest for biological applications (magnetic hyperthermia) as well as for heterogeneous catalysis. This study shows why FeC NPs of similar structures and static magnetic properties display radically different heating power (SAR from 0 to 2 kW.g -1 ). By combining results from Transmission Electron Microscopy (TEM), Dynamic Light Scattering (DLS) and static and time-dependent high-frequency magnetic measurements, we propose a model describing the heating mechanism in FeC nanoparticles. Using, for the first time, time-dependent highfrequency hysteresis loop measurements, it is shown that in the samples displaying the larger heating powers, the hysteresis is strongly time dependent. More precisely, the hysteresis area increases by a factor 10 on a timescale of a few tens of seconds. This effect is directly related to the ability of the nanoparticles to form chains under magnetic excitation, which depends on the presence or not of strong dipolar couplings. These differences are due to different ligand concentrations on the surface of the particles. As a result, this study allows the design of a scalable synthesis of nanomaterials displaying a controllable and reproducible SAR.

INTRODUCTION

 As a consequence, variations in the reaction conditions may lead to variations in the heating power of the MNPs, which emphasizes the need for robust synthetic methods for MNPs production.

 The

magnitude of magnetic interactions (dipolar couplings) between the NPs is crucial to allow this phenomenon taking place, since if the NPs are too far away to interact and to be organized in the presence of the magnetic field, their heating power remains low. Thus, the heating power increases with the concentration of NPs until reaching an optimal point, after which higher concentrations favor the presence of strong demagnetizing interactions resulting into a drop of the heating power. 20 One additional consequence of dipolar interactions is that, at a given concentration, an agglomeration of MNPs reduces their heating due to an absence of mobility 25 and/or to the fact that the overall magnetic moment of agglomerates is reduced due to their magnetic flux closure configurations, reducing their tendency to interact with their neighbors. 26,27 Thus, in samples without agglomerates and displaying large heating powers, the heating power is expected to increase with time following the dynamics of the NPs self-organization upon application of an alternating magnetic field. On the contrary, in samples with agglomerates, heating power is expected to remain small and constant with time. However, to the best of our knowledge, direct experimental proofs of these expected effects have never been reported. This will be done in the present work and monitored by a technique never used in this context: time-dependent highfrequency measurements.

Different MNPs have been synthesized for applications in magnetic hyperthermia. Among them, the most widely used are iron oxide-based systems due to their good biocompatibility. 26,[28][29][30][31][32][33] However, for magnetically-induced catalysis where higher temperatures are required, iron carbide NPs (FeC NPs) can be better candidates. especially since these materials are ferromagnetic at room temperature and display moderate coercive fields and good colloidal stability. 2 Different preparation methods of FeCx NPs for applications in catalysis have been reported. [34][35][36][37][38][39] Over the past years, our group has developed the synthesis of monodisperse Fe(0) NPs, either pure or bimetallic [40][41][42][43] and their transformation into monodisperse iron carbide NPs. 5,34,44 High values of SAR of ca. 750 Wˑg -1 (100 kHz, 47 mT) were achieved for Fe(0) NPs, 45 but their transformation into FeC NPs led to much enhanced SAR as well as better catalytic activity and air stability. Indeed, we recently reported the synthesis of FeCx NPs displaying heating properties in organic solvents that surpass those of any previously described material. 5,44 The NPs were synthesized by carbidization of pre-formed monodisperse Fe(0) NPs using a mixture of syngas (CO/H2). We demonstrated that their very high SAR was related to the presence of the crystalline Fe2.2C phase. 5 However, different values of SAR could be obtained upon slightly modifying the reaction and purification conditions.

In this paper, we demonstrate that standard static magnetic measurements cannot discriminate between heating and non-heating Fe2.2C NPs. In our case, the heating properties can only be understood in terms of the presence or not of strong dipolar couplings resulting from the organization or agglomeration of the nanoparticles, which in turn depend on the surface ligands concentration. Time-dependent high-frequency magnetic measurements were employed for the first time to characterize the behavior of NPs assemblies under alternating magnetic fields. In association with dynamic light scattering they were found to be the key techniques to unveil the origin of the differences in the heating properties of the samples. These findings constitute a novel insight into the relationship between the synthesis of NPs and their final heating properties and allowed to design robust procedures for the production of FeC NPs based nanomaterials displaying high and reproducible heating powers.

RESULTS AND DISCUSSION

Heating and non-heating Fe2.2C NPs obtained by carbidization of Fe(0) NPs.

In our previous work, we explored the synthesis of Fe2.2C NPs through carbidization of preformed monodisperse Fe(0) NPs under a syngas mixture in mesitylene. 5 In that preliminary study, the carbidization was performed starting from Fe(0) NPs of 12.5 nm. The Fe(0) NPs had been isolated by post-synthesis decantation and washing, after which a powder with a Fe content of ca. 50 wt% was obtained and carbidized during 96 h to give Fe2.2C (FeC-ref) NPs with unprecedent SAR (see Table 1). The excellent heating power of the Fe2.2C NPs obtained after carbidization incited us to optimize and scale-up the reaction. Thus, in order to purify more extensively the surface of the particles, a more efficient washing of the powder led to Fe(0) NPs Fe-1, with a Fe content of 78% determined by thermogravimetric analysis (TGA). Surprisingly, even though the carbidization of these NPs was performed under similar conditions to those described in our previous work (reaction time of 96 hours), a strong agglomeration of the FeC NPs (NPs FeC-1 in Table 1, see also Figure 1-a and S1) was observed by transmission electron microscopy (TEM).

The SAR of FeC-1 measured by calorimetry experiments had a value of 0 W•g -1 when applying an alternating magnetic field of µ0Hrms of 47 mT at a fixed frequency of 93 kHz, (Hmax = √2Hrms) contrasting with the excellent SAR found under the same conditions for NPs FeC-ref prepared by the same methodology (ca. 2000 W•g -1 , see Table 1). In addition, NPs FeC-1 were not dispersible in mesitylene, in contrast to the NPs FeC-ref prepared in our previous work. Furthermore, the NPs rapidly precipitated even after external addition of solubilizing ligands and sonication at room temperature. Finally, the size of the FeC NPs FeC-1 after carbidization was only marginally increased compared to the starting Fe(0) NPs (from 12.5 ± 1.3 nm to 13.7 ± 1.6 nm). Extending the carbidization time to 7 or 8 days induced the coalescence of the NPs to give large nanoobjects of more than 200 nm observed by TEM. These results suggest a noninnocent role of the ligands in the carbidization process.

Analysis of the Mössbauer spectrum (see Fig S2 ) showed the presence of the heating carbide phase Fe2.2C in FeC-1 with a concentration of 59.6%, together with some Fe5C2 (28.1%) and Fe(0) (8.3%). Although the relative amount of the heating phase Fe2.2C was lower than in FeCref NPs (71.2% of Fe2.2C in the latter), the absence of any heating capacities in FeC-1 cannot be justified by this difference. Powder X-Ray diffraction (XRD) analysis of FeC-1 only showed the peaks related to the Fe2.2C phase (see Fig. S3). Similarly, the saturation magnetization (Ms, 170 A•m 2 •kg -1 ), remnant magnetization (Mr, 23 A•m 2 •kg -1 ), coercive field (Hc, 0.051 T) and High Resolution TEM (HRTEM) and High-Angle Annular Dark-Field Scanning Transmission Electron microscopy (HAADF-STEM) analyses of the non-heating FeC-1 NPs showed that the agglomerates formed by several NPs are surrounded by a thin amorphous layer composed of a light element (see Figure S5). Energy-dispersive X-ray EDX spectroscopy also evidenced the presence of a shell containing silicon in addition to oxygen (see Figure 2b-f). This amorphous silicon-oxygen embedding layer is likely to be amorphous silica or a polysiloxane derivative, the presence of which could be related to the lack of solubility of NPs FeC-1 in mesitylene and their assembly into the large aggregates as it would reduce the interactions between the agglomerates and the ligands added to re-disperse the NPs. XPS spectroscopy confirmed the presence of Si at binding energies of 102 eV for the Si 2p3/2 orbital, which are characteristic of silicates or organosiloxanes (see Fig S6). Concerning its origin, the only possible source of Si in the carbidization process is hexamethyldisilazane (hmds) resulting from the reaction with hydrogen of the corresponding amido ligand present in {Fe[N(SiMe3)2]2}2, the organometallic precursor of the starting Fe NPs (see experimental section for more details). When carbidization progresses, water molecules are released into the solution as a reaction product 44 and at low concentration of other ligands they can react with the Si moiety, leading to the formation of the amorphous embedding layer.

The role of this layer in the on-off behaviour of the heating power has then to be clarified. It has been previously shown that the presence of more than 2 wt% of SiO2 deteriorates the magnetic permeability of a hard magnetic material such as bulk barium ferrite, 46 but here, this is not the case since the Ms of the FeC-1 NPs determined by VSM is very similar to that measured for the FeC-ref NPs (see Table 1). Then, the on-off behaviour of the heating power has to be related to the NPs agglomeration, in which the silicon containing layer can also play a role by preventing chain formation. There are two convergent explanations to this fact: (i) when the NPs are aggregated, the proximity between the magnetic domains induces strong demagnetizing interactions between them, typically dipole-dipole interactions 18,20,47 and (ii), when the FeC NPs are agglomerated, their organization in the magnetic field into chains is not possible. Such chains are necessary to be formed during the hyperthermia experiments in order to induce the presence of a uniaxial anisotropy that enhances the NP heating power. 19,22 For this phenomenon to happen, a relative free mobility of the NPs is required.

We noticed that after 96 hours of carbidization, 2 bars of syngas were consumed for the preparation of the FeC-1 NPs, whereas for the synthesis of the heating FeC-ref NPs the pressure drop observed at the end of the reaction was of only 0.8 bars. We assume that due to the presence of less organic ligands in the Fe-1 NPs herein used, an easier surface access may accelerate CO hydrogenation, leading to more syngas consumption and justifying the higher pressure drop. When carbidization of Fe-1 NPs was allowed to proceed only for 48 hours, the syngas consumption was of 0. their SAR measured at µ0Hrms of 47 mT with a f of 93 kHz was 1260 W.g -1 (see Table 1 and Interestingly, despite their better heating properties, the degree of carbidization of the FeC-2

NPs was lower than for FeC-1, as deduced from the lower content of the Fe2.2C and Fe5C2 phases in the Mössbauer spectroscopy (see Table 1). This fact also indicates that the absence of heating power in NPs FeC-1 is not related to the carbidization degree and raises the question of the heating mechanism. From the above data we can conclude that the low ligand content of the Fe-1 NPs modifies the kinetics of the carbidization process and that the FeC-1 NPs present different surface properties and tend to agglomerate. Indeed, XPS spectrometry revealed that NPs FeC-2 were free of Si, contrary to FeC-1. XPS also showed that NPs FeC-2 were richer

Magnetic properties c

Ms

= 170 A•m 2 •kg -1 Mr = 23 A•m 2 •kg -1 Hc = 0.051 T χm = 0.043 Ms = 168 A•m 2 •kg -1 Mr = 28 A•m 2 •kg -1 Hc = 0.055 T χm = 0.043 Ms = 170 A•m 2 •kg -1 Mr ~ 30 A•m 2 •kg -1 Hc = 0.103 T χm ~ 0.04 Ms = 163 A•m 2 •kg -1 Mr ~ 25 A•m 2 •kg -1 Hc = 0.050 T χm ~ 0.04 Mössbauer spectroscopy d Fe2.2C (60 %) Fe5C2 (28%) Fe(0) (8%) Fe2.2C (52%) Fe5C2 (32%) Fe(0) (12%) Fe2.2C (71%) Fe5C2 (24%) Fe(0) (4%) Fe2.2C (55%) Fe5C2 (23%) Fe(0) (18%)
in Fe(0) at their surface and less oxidized than NPs FeC-1, which may be again linked to the presence of the Si containing layer in the latter (see Fig S6).

In order to re-disperse the agglomerates, a sample of 50 mg of the FeC-1 NPs was heated at 150 °C in mesitylene in the presence of 0.2 equivalents of palmitic acid for 2 hours. A darkbrown suspension that contained the FeC-3 NPs was recovered after the reaction. Although the suspension was stable for 1 hour and only thereafter the NPs started to precipitate, large agglomerates were still observable by TEM (see Figure S12) together with few isolated NPs. and The FeC-3 NPs displayed a very low SAR value of ca. 20 W.g -1 in hyperthermia experiments at µ0Hrms = 47 mT and f = 93 kHz. The FeC-3 NPs were then characterized by VSM and XRD without showing significant differences from FeC-1 and FeC-2 NPs (see Fig. S13-15). Moreover, the hysteresis cycles for NPs FeC-3 and FeC-2 are identical. Then, the FeC-3 NPs were further characterised by HAADF-STEM and EDX. No significant differences in the surface composition of the FeC-3 NPs and the heating FeC-2 NPs were found by EDX mapping. However, upon observing TEM micrographs at low magnification, it is difficult to determine if the agglomeration degree of the two samples FeC-2 and FeC-3 is significantly different (see Figure S14).

Understanding the heating mechanism of FeC NPs in solution.

Thanh and co-workers have observed a correlation between the heating properties of multi-core iron oxide NPs and their dispersion in solution. 27 The authors concluded that a better magnetic heating power was associated to less inter-particle interactions as deduced by Dynamic Light Scattering (DLS). Thus, DLS measurements in mesitylene solution evidence for FeC-3 the presence of large agglomerates displaying a DH of 590 ± 380 nm (38% population) and 121 ± 40 nm (62% population) (see Figure 2). In contrast, for FeC-2 displaying a good heating power, a single distribution was observed in the DLS with a DH of 14.3 ± 1.0 nm. This indicates that inter-particle interactions in solution are a determining factor for explaining the heating properties as previously proposed. In order to confirm this hypothesis, we studied the intrinsic magnetic properties of FeC-2 and FeC-3 by dispersing the nanoparticles in tetracosane at 60 °C and performing magnetic measurements. Tetracosane solidifies at low temperature and freezes the position of the particles in the medium. This then limits the magnetic interactions between the particles (see Table 2). The FeC-1 NPs could not be studied as they were not dispersible due to the presence of the larger agglomerates. A zero-field cooled -field cooled (ZFC/FC) experiment with a µ0H0 of 10 mT was run, together with the hysteresis loops measurements at 300 K and 5 K (see Fig S16-18). Even though the NPs were diluted in tetracosane to very low Fe concentrations (Fe wt% of ~ 0.01%), in both cases the NPs were far from being well-dispersed and strong interactions between the particles still remained. Thus, the apparent superparamagnetic blocking temperature (Tb) at which the ZFC and FC graphs matched was above 300 K, which we attribute to the agglomerates present at different degree in each sample. 48 Higher values of Hc and Mr, as well as a saturation of the loops at lower amplitudes of magnetic field were found for FeC-2 in the hysteresis loops at 5 K after dilution in tetracosane, which agrees with weaker interactions in this sample. ZCF/FC curves indicate that both samples are aggregated but no qualitative conclusion on the difference between the two samples can be provided.

To gain a more precise understanding of the magnetic properties of this system, we turned to measurements under dynamic conditions which are more representative of the situation under hyperthermia conditions. We measured the high-frequency hysteresis loops of both samples under an alternating magnetic field of µ0Hrms = 37 mT and 50 kHz using a previously described set-up. 49 The high-frequency hysteresis cycles for the FeC-2 and FeC-3 were completely different (see Fig S19 ): FeC-2 displayed a fully opened hysteresis loop with a maximum magnetization at 37 mT of 84 A•m 2 •kg -1 , whereas FeC-3 displayed a very small response to the magnetic excitation, with an almost reversible response and a maximum magnetization at 37 mT of 28 A•m 2 •kg -1 . The hysteresis area is 23 times larger for FeC-2 than for FeC-3, in agreement with its much more efficient heating. These results demonstrate that although the magnetic properties of the NPs under static conditions are similar for both systems, their dynamic magnetic properties are not comparable. a Measured at 47 mT and 93 kHz. b Ms: Saturation magn. Hc: Coercive field. Measured by VSM at 3T and 300 K. c Ms: Saturation magn. Hc: Coercive field. Measured by VSM at 3 T and 300 K after dilution in tetracosane (Fe wt% of ~ 0.01%). d Mmax: Magnetization at maximum magnetic field (µ0Hrms of 37 mT). Hc: Coercive field. Cycles were measured in solution at 300 K (10 mg of Fe2.2C NPs in 0.5 mL of Mesitylene).

To have a deeper insight into the origin of the difference in heating between the two samples, high-frequency hysteresis loops for FeC-2 and FeC-3 were recorded as a function of time using the same set-up. A very different behaviour was found between the two samples (see Fig 3).

Whereas the hysteresis loops for NPs FeC-2 open with time, reaching their maximum value of area after 60 seconds, the hysteresis cycles of the non-heating NPs FeC-3 remained closed. Fig. S20 displays the evolution of the hysteresis area as a function of time for the two samples.

Interestingly, the initial heating power of both samples is quite small and comparable, since the hysteresis area is approximately 0.41 and 0.14 Jˑkg -1 for FeC-2 and FeC-3 respectively. However, in the case of FeC-2, this value increases up to 7.23 Jˑkg -1 after 60 s without significant evolution afterwards, whereas it remains in the same range for FeC-3. Similarly, Fig

S21 displays the values of Mmax and normalized χm as a function of time for NPs FeC-2 and

FeC-3. Whereas the Mmax and χm values remain almost constant for NPs FeC-3, Mmax values

get increased by a factor of 7 and χm values by a factor of 10 for NPs FeC-2 after 60s of exposure to the alternating magnetic field. This is the clear sign that, in sample FeC-2, the NPs organize into chains which can be macroscopically observed with a naked eye. This organization induces an effective uniaxial anisotropy along the magnetic field direction, which enhances the heating power of the material. In addition, the area of the high-frequency hysteresis loop after 10 s of exposure has a value of 4.04 Jˑkg -1 (see also Fig S20 ), which is in of the same order of magnitude as the value measured by calorimetry after the same time of exposure at µ0Hrms of 33 mT and a f of 93 kHz (5.67 Jˑkg -1 ). It should be noted that, for experimental reasons, the concentration in mesitylene of NPs for SAR measurements by calorimetry experiments is ca. 20 mgˑmL -1 , whereas in the dynamic hysteresis measurements it is ca. 10 mgˑmL -1 . Qualitatively, we observed the rapid (1-2 s) formation of chains when the experiments were performed at 5 different concentrations of NPs FeC-2 (from 2 to 40 mgˑmL -1 ). The discrepancy between the measurements for FeC-2 and FeC-3 might be due to the fact that larger aggregates are not able to modify their structures to form chains, or that, due to their flux closure magnetic states and the larger inter-aggregate distance, they interact less with each other. To the best of our knowledge, this is the first time that high-frequency hysteresis loops as a function of time are measured for characterizing the heating power of MNPs. First, it clearly evidences that heating power can be a time-dependant property, which was overlooked in previous works. Second, the results obtained on these two samples illustrate that timedependent high-frequency hysteresis loops constitute a precious tool to understand the origin of heating power in solution of MNPs: whereas the two samples we had measured displayed similar magnetic properties under VSM, their heating properties were very different, which is related to the possibility for one of them to organize in solution.

Finally, we carried out an additional experiment to observe the formation of chains under hyperthermia conditions. When a drop of the above-mentioned solutions was deposed on a grid in the presence of an alternating magnetic field of µ0Hrms = 33 mT and 93 kHz and left to evaporate to dryness, the formation of chains of size within the range of 0.5-5 µm that were aligned with the magnetic field direction was clearly observable in the TEM micrographs of We conclude this part by presenting a summary of our view on the mechanism at the origin of heating in these samples. The overall hysteresis loop of assemblies of MNPs is the combined result of i) their physical rotation under the influence of the torque or of thermal agitation, and ii) their magnetization reversal. In our case, we consider that agglomerates of MNPs do not rotate under the influence of the magnetic field because they are rather large objects with a very weak remnant magnetization so the torque they undergone is extremely small. In the case of chains, however, we cannot exclude that part of the reversal is due to their physical rotations, but there is so far no suitable model permitting to easily calculate the reversal time of an anisotropic assembly of nanoparticles. Moreover, even if this model existed, it would require knowing precisely the length and width of the chains during the application of the magnetic field. In any case, whatever the underlying mechanism, Fig. 3(a) illustrates well that the coercive field of the chains is close but below the applied alternating magnetic field, which explains why FeC particles display large heating power upon magnetic moment reversal.

Ligand effects on the heating power of Fe2.2C NPs: reproducible and scalable preparation of nanoparticles displaying high heating power

The previous results demonstrate that the initial ligand concentration seems to play a crucial role in the carbidization process. This concentration depends on the washing process of the nanoparticles. In this work, we have tried to remove as much ligands as possible before reacting the particles with syngas. This led to Fe-1 NPs with a Fe concentration of 70-80 wt%. Starting from this sample, and in order to slow down the carbidization process in a controlled way avoiding at the same time the formation of the agglomerates observed in NPs FeC-1, a known amount of ligands was added to the NP suspension in mesitylene. In this way, the relation between the nature of the ligands and their concentration in the reaction medium with the SAR of the FeC NPs obtained could be studied.

Four different carbidization experiments were carried out (see experimental section for details), starting from the NPs Fe-1 with a Fe content of 73 wt% determined by TGA. In these experiments, various amounts of a 1:1 (by weight) mixture of hexadecylamine (HDA) and palmitic acid (PA), the ligands used for the synthesis of the Fe(0) NPs, were introduced into the reactor before reacting with syngas. The resulting NPs FeC-4 to FeC-7 were thoroughly characterized and their characteristics are summarized in Table 3 (see also Figures S23-28).

When carrying out the reaction after addition of solubilizing ligands, all FeC NPs samples could be dispersed in mesytilene and no large agglomerate was observed on the TEM images. In all cases the resulting NPs displayed heating properties in mesitylene solution and a correlation between the SAR and the amount of ligands added before the carbidization process was established. The SAR of the NPs progressively increased when increasing the amount of the 1:1 HDA/PA mixture (see Figure 5). A maximum value of 1700 W•g -1 at a field of µ0Hrms of 47 mT with a frequency of 93 kHz was obtained for NPs FeC-6 when the reaction was performed starting from a solid sample displaying a 37 wt% Fe content. Further lowering the Fe content down to 24 wt% led to a decrease of the SAR value of NPs FeC-7 to 670 W•g. -1 It is also interesting to note that syngas consumption at the end of the reaction progressively decreased from FeC-4 to FeC-7, meaning that CO hydrogenation is slowed down when the amount of ligands present in the reaction medium was increased. This observation can be related to the surface accessibility of the nanoparticles. In addition, the crystallite sizes of the Fe2.2C phase determined by the Scherrer equation from the XRD were slightly larger upon increasing the amount of ligands. In all cases, extending the carbidization time to 6 or 7 days did not lead to any notable difference in the heating power of the FeC NPs.

Then, two carbidization experiments were performed under the reaction conditions in which the highest SAR was obtained, i.e. starting from a Fe content of 37%, but adding only PA (FeC-8) or HDA (FeC-9) as ligands. NPs FeC-8 and FeC-9 were characterized by all the above-mentioned techniques (see Fig. S29-33). The nature of the ligand was shown to have a strong influence on the CO hydrogenation rate as deduced from the pressure drop in both experiments at the end of reaction. More precisely, the pressure drop observed was 0.7 bar for FeC-8 and of >1.5 bar for FeC-9, which can be related to the stronger coordination of PA to the NPs as compared to the HDA. In addition, PA was found to produce FeC NPs displaying higher heating power. Thus, the SAR measured for the FeC-8 NPs was 2020 W•g -1 whereas it was only 610 W•g -1 for the FeC-9 NPs when applying a field of µ0Hrms of 47 mT and a frequency of 93 kHz (see Figure 5 and Experimental section).

The agglomeration degree was again found to play a role on the heating power, since the FeC-9 NPs were much more agglomerated than the FeC-8 NPs in agreement with the much better coordination properties of PA compared to HDA. In addition, on the XRD diffractogram of the FeC-9 NPs, some peaks corresponding to iron oxide were observed. This can be attributed to the competition between water resulting from the carbidization process and the amine ligand for iron coordination, competition which is not present when PA is used. NPs FeC-8, were also studied by Mossbauer spectroscopy (see Fig S31). The percentages of the Fe2.2C (69.6%) and Fe5C2 (24.8%) phases were now comparable to those obtained by Bordet et al. under similar conditions (71.2% of Fe2.2C and 24.4% of Fe5C2), which displayed also a comparable value of SAR.

Finally, high-frequency hysteresis loops as a function of time were acquired for NPs FeC-5 and FeC-8 (see Fig S34 -37) at µ0Hrms of 33 mT and a f of 50 kHz at a concentration of ca. 10 mgˑmL -1 . The area of the cycles for NPs FeC-5 after 40-60 s were comparable to those found for NPs FeC-2 (ca. 6-7 mJ•g -1 ), which agrees with their similar values of SAR obtained by calorimetry for these samples. The dynamic of chain formation was also comparable, and the hysteresis loops began to open after 1-2 s of exposure to the alternating magnetic field.

However, for the NPs FeC-8 displaying the highest value of SAR, the area of the hysteresis loops reached a maximum of 14 mJ•g -1 after 5-10 s of field exposure, but rapidly decreased because the chains were coming out of the solution as observed in the calorimetry experiments and thus their heating power could not be measured anymore (see Fig S38). Furthermore, chain formation was faster in this case, and the hysteresis loops began to open very rapidly, 0.1 s of exposure to the magnetic field. The mobility in solution of well-dispersed FeC NPs is an important factor to obtain high SAR values. The nature of the ligand also influences the preparation of FeC NPs. PA coordinates more strongly to the surface of the Fe(0) NPs than HDA. 42 The oxidation of the NPs FeC-9

observed in the XRD is in agreement with this hypothesis. If HDA gives weaker interactions with the NPs, they are more prone to react with the water that is generated in the carbidization process leading to partial oxidation of the surface. In addition, PA may assist the shift of the surface-coordinated hmds by protonation, hence helping to avoid the formation of the Si containing amorphous layer that embeds the agglomerates composing the non-heating NPs FeC-1. Understanding the influence of the reaction conditions allowed preparing reproductively and in a scalable way NPs FeC-8, which display the highest heating power in this work (ca. 2000 Wˑg -1 ) at a field of µ0Hrms of 47 mT with a f of 93 kHz.

Scheme 1. Schematic representation of the role of the ligands into the heating properties of the Fe2.2C NPs.

CONCLUSIONS

In this work, in order to scale-up the reproducible synthesis of highly heating magnetic nanoparticles and to understand their heating mechanism, we have extended our previous studies concerning the synthesis under syngas of FeC NPs by carbidization of pre-formed Fe(0)

NPs. We have studied their structure and their heating properties in detail. We have discovered a surprising effect of the extended washing of the starting Fe(0) nanoparticles. When the nanoparticles surface is depleted of its ligands after extended washing, the particles are strongly and irreversibly agglomerated and cannot align in an alternating magnetic field. The result is the absence of heating properties of these nanomaterials. Control over the ligands concentration on the surface of the NPs is therefore a pre-requisite to prepare nanoparticles displaying a controllable and reproducible heating power. In parallel, this work explores the heating mechanism in magnetic nanoparticles, demonstrates the importance of dynamic magnetic measurements and brings forward a new indication of the importance of dipolar interactions, here specifically in solution, for controlling the heating power of magnetic materials.

Finally, when controlling the Fe concentration at the beginning of the reaction, the carbidization process showed good reproducibility and scalability. We believe that this work represents an important advance in the synthetic control of well-dispersed homogeneous magnetic NPs, here FeC NPs, displaying reliable heating properties in the objective of future applications.

CONCFLICTS OF INTEREST

The

  magnetic susceptibility (χm, 0.043) of the FeC-1 NPs measured by vibrating sample magnetometry (VSM, see Fig S4) did not reveal any significant differences in the Ms, Mr or χm from the values found in the FeC-ref NPs previously reported by Bordet et al., although a higher value of Hc was obtained in that case (Ms = 170 A•m2•kg-1, Mr ~ 30 A•m2•kg-1, Hc = 0.103 T and χm ~ 0.04 at 300 K). All these results indicate that these structural differences are not enough to explain the lack of heating properties in FeC-1 NPs.

Figure 1 .

 1 Figure 1. TEM micrographs of (a) the agglomerated NPs FeC-1 after carbidization of 96 h and

  8 bars, and the FeC-2 NPs obtained were dispersible in mesitylene as evidenced by TEM micrographs (see Fig 1b and Fig S7). In contrast to FeC-1,

Fig S8). The FeC- 2

 2 Fig S8). The FeC-2 NPs were also characterized by the same techniques as the FeC-1 NPs (see

Table 1 . 5 a

 15 Comparison between the NPs FeC-1 and FeC-2 prepared in this work by variation of the reaction time with the NPs FeC-ref prepared by Bordet et al. Deterimined by TGA. b Measured at 47 mT and 93 kHz. c Determined by VSM and measured at 300 K. d Measured at 4 K.

Figure 2 .

 2 Figure 2. Size distribution of the hydrodynamic diameters measured by DLS for the NPs FeC-

Figure 3 .

 3 Figure 3. High-frequency (50 kHz) hysteresis loops as a function of time for the heating NPs

FeC- 2

 2 (see Fig 4. and Fig S22). However, in the case of the non-heating NPs FeC-3, no anisotropic orientation along the magnetic field could be observed and the NPs remained in the form of big amorphous agglomerates.

Figure 4 .

 4 Figure 4. TEM micrographs of the NPs (a) FeC-2 and (b) FeC-3 when a drop of the

Figure 5 .

 5 Figure 5. (a) Graphical representation of the SAR measured at a field of µ0Hrms of 47 mT with

  authors declare no conflict of interest AKNOWLEDGMENTS The authors thank ERC Advanced Grant (MONACAT 2015-694159) for financial support, J.-F. Meunier for Mössbauer measurements and Jérôme Esvan and Simon Tricard for the XPS measurements. . Ortgies, F. J. Teran, U. Rocha, L. de la Cueva, G. Salas, D. Cabrera, A. S. Vanetsev, M. Rähn, V. Sammelselg, Y. V. Orlovskii and D. Jaque, Advanced Functional Materials, 1704434-n/a.

Table 3 .

 3 Comparison between the heating NPs FeC-2 and non-heating NPs FeC-3.

	Properties	NPs FeC-2	NPs FeC-3
	Size (nm)	14.4 ± 1.3	12.1 ± 1.1
	SAR (W•g -1 ) a	1260	20
	Ms: (A•m 2 •kg -1 )	Ms: 170; Hc:	Ms: 168; Hc:
	Hc: (T) b	0.051	0.055
	Ms: (A•m 2 •kg -1 ) Hc: (T) in tetracosane c	Ms ~ 151; Hc: 0.020	Ms ~ 135; Hc: 0.025
	Mmax: (A•m 2 •kg -1 )
	Hc: (T) at µ0Hrms of	Mmax ~ 84; Hc:	Mmax ~ 28; Hc:
	37 mT and f of 50	0.037	0.005
	kHz d
			590 ± 380
	DH by DLS (nm)	14.3 ± 1.0	(38% in number) 121 ± 40
			(62% in number)

Table 3 .

 3 Comparison between the characterizations of NPs FeC-4 to FeC-9 prepared by carbidization of the NPs Fe-1 after addition of different amounts of a 1:1 mixture in mass of HDA/PA, PA or HDA.

	FeC NPs	FeC-4	FeC-5	FeC-6	FeC-7	FeC-8	FeC-9
	Wt% of Fe in starting solid a	52%	46%	37%	25%	37% (PA)	37% (HDA)
	ΔP after reaction (bar)	1.6	1.2	0.8	0.5	0.7	2.0
	Size (nm)	13.7 ± 1.4 nm	13.8 ± 1.4 nm	14.9 ± 1.5 nm	12.6 ± 1.5 nm	14.6 ± 1.3 nm	13.7 ± 2.0 nm
	XRD crystallite Size (nm)	9.3	9.9	10.2	10.7	9.9	8.6 (Oxidized)
	SAR (W•g -1 ) b	970	1130	1700	670	2020	610
	Mössbauer spectroscopy	-.	Fe2.2C (56%) Fe(0) (5%) Fe5C2 (29%)	-	-	Fe2.2C (70%) Fe(0) (2%) Fe5C2 (25%)	-

a After addition of a 1:1 PA/HDA mixture in FeC-4 to FeC-7, PA in FeC-8 or HAD in FeC-9. b Measured at µ0Hrms of 47 mT and a f of 93 kHz.