Cyril Roscian

Florian Praden

Jean-Max Dutertre

Jacques Fournier

Assia Tria

SECURITY CHARACTERISATION OF A HARDENED AES CRYPTOSYSTEM USING A LASER

Keywords: Advanced Encryption Standard (AES), fault attacks, laser, security characterization, DFA

The AES is a standard encryption algorithm used in numerous cryptographic systems like smart cards, TPMs as well as in protocols like WPA2 or OpenSSL. Measuring the robustness of AES implementations against physical attacks is of utmost importance in order to guarantee the security of the system into which the AES is used. In this article, we describe how a hardware AES, embedding countermeasures against physical attacks, has been characterized using a laser. With the latter, we tried to implement a class of physical attacks called fault attacks which, when successful, could have helped us retrieve the secret key used by the AES module. Our experiments have allowed us to validate the efficiency of some of the countermeasures implemented in this AES implementation and have given us hints on how to further improve them.

Introduction

Security, through the authenticity, confidentiality and integrity of communication systems, has become an essential component of all electronic systems. The vulnerability to attacks of the devices implementing the cryptographic algorithms (such as smart cards) has become a critical issue. Some malicious means or "physical attacks" could be used to retrieve sensitive information such as encryption keys and therefore lower the security of the whole protected transmission chain of information. There are three types of physical attacks: invasive attacks, which cover all the techniques based on the "physical" analysis, modification [START_REF] Anderson R | Low Cost Attacks on Tamper Resistant Devices[END_REF] and probing [START_REF] Handschuh | Probing Attacks on Tamper-Resistant Devices[END_REF], SCHMIDT 2009[START_REF] Gammel | On the Duality of Probing and Fault Attacks[END_REF] of integrated circuits (IC) by an invasive method [START_REF] Kömmerling | Design principles for tamper-resistant smartcard processors[END_REF], KOEUNE 2005); observation or passive attacks which exploit the fact that some physical characteristics such as power consumption [START_REF] Kocher | Differential Power Analysis[END_REF][START_REF] Lu | Principles on the security of AES against first and second-order differential power analysis[END_REF] , electromagnetic radiation or the duration of computation depends on the chip's internal calculations (KOEUNE 2005); perturbation or fault attacks, which are based on changing environmental conditions of the chip. The latter ones are the most complex to implement as various and complex parameters must be taken into account such as the timing (i.e. the synchronization between the injection time and the calculation), the localization and the power level.

In this paper we mainly discuss about characterisations based on fault attacks. One of the goals of such attacks could be to reveal the secret keys of a cryptographic device based on techniques like Differential Fault Analysis (DFA) [START_REF] Biham | Differential Fault Analysis of Secret Key Cryptosystems[END_REF][START_REF] Boneh | On the Importance of Checking Cryptographic Protocols for Faults[END_REF]. Several methods can be used to induce faults into cryptographic integrated circuits: use of a laser [START_REF] Skorobogatov | Semi-Invasive Attacks --A New Approach to Hardware Security Analysis[END_REF], BAR-EL 2006), voltage pulses [START_REF] Blömer | Fault Based Cryptanalysis of the Advanced Encryption Standard (AES)[END_REF], clock glitches (AGOYAN 2010) or electromagnetic (EM) disturbances [START_REF] Schmidt | Optical and EM Fault-Attacks on CRT-based RSA: Concrete Results[END_REF].

The use of lasers is one of the most effective techniques: lasers allow a good reproducibility, an accurate control on the timing of the injection (the instant of firing and the duration of the pulse) and a precise focalization (the location and the penetration of the fault). Consequently, lasers generate precise local effects onto the IC, thus 'protecting' the remaining chip against unwanted perturbations. Several theoretical attacks against cryptographic algorithms are based on such models of fault injection methods [START_REF] Piret | A Differential Fault Attack Technique Against SPN Structures, with Application to the AES and KHAZAD[END_REF][START_REF]DFA on AES[END_REF], 2010, MORADI 2010). Although protections against these kinds of attacks exist [START_REF] Yen | Simple Error Detection Methods for Hardware Implementation of Advanced Encryption Standard[END_REF], advanced methods combine semi-invasive attacks and power or EM analysis [START_REF] Moradi | On the power of fault sensitivity analysis and collision side-channel attacks in a combined setting[END_REF].

We hereby describe the security assessment of an AES hardware chip done to validate the efficiency of the embedded countermeasures that could be incorporated into the AES module of the SECRICOM's secure docking station (SDM). The SDM contains authentication keys and access rights associated with each user of the SECRICOM secure communication network. To access such keys, a secure channel, based on AES encryption, has to be established between the SDM and the host. Hence, guaranteeing the resistance of the AES used against physical attacks helps in hardening the security chain of SECRICOM's communication infrastructure.

The outline of this paper is as follows. First, we recall the physical phenomena occurring when a laser is used to inject a fault into an IC. Then we describe the AES implementation used as target of our laser-based characterization. Third, we provide a description of the laser-based tests made on the targeted AES. With attacks like those described in [START_REF] Piret | A Differential Fault Attack Technique Against SPN Structures, with Application to the AES and KHAZAD[END_REF][START_REF]DFA on AES[END_REF] in mind, we tested the countermeasure in a "black box" approach (i.e. without using the knowledge we had of the implemented countermeasures). We also tested the AES in a "white box" approach where we tried to take advantage of the knowledge we had of the implemented countermeasures in order to circumvent them. Finally we provide a discussion about the results that we obtained and the conclusions that could be drawn from them.

Physical phenomena induced by lasers into an IC

When an IC is exposed to a laser beam, two phenomena have to be considered: the photoelectric effect appearing and the fact that some parts of the IC are more sensitive than others. When a laser beam strikes the silicon and that the photon's energy is higher than the Silicon's band gap, electron-hole pairs are created.

In general, these pairs recombine and there is no effect on the IC but under specific conditions, some undesired effects can appear. We shall focus on one of those effects called the Single Event Effect (SEE).

Single Event Effect

A Single Event Effect can appear when the electron-hole pairs created by the laser beam are drifted in opposite directions by the electrical field instead of immediately recombining. The consequence of that is the creation of a transient current as illustrated in Figure 1. (hence the extension of the electric field) along the laser beam and within a few picoseconds, charges are collected giving a current peak. In a second time, the rest of the charges are collected in a longer diffusing scheme. Figure 2 shows the transient current associated with the two phenomena of collection as given in [START_REF] Wang | Single Event Upset: An Embedded Tutorial[END_REF].

Sensitive zones

In CMOS technology, some parts are more sensitive than others to SEEs. Indeed, to create an SEE, and then a fault, a strong electric field is needed. The reverse biased PN-junctions of the chip provide this required electric field. The position of these junctions can change, depending on the value of the data manipulated.

A good example to illustrate this data dependency is the CMOS inverter. The first case is a high state on the inverter's input: the NMOS transistor is in ON state, its drain is grounded, the source and the bulk too and there is no reverse biased PN-junction. The PMOS transistor is in OFF state, then its source and the N well are in high potential, but its drain is grounded giving rise to a reverse biased junction. Hence, the drain of the PMOS transistor becomes sensitive to a laser shoot. In the same way, with a low state on the inverter's input, the drain of the NMOS transistor becomes sensitive to a laser shoot.

The hardened AES test chip

In this section, we describe the chip used as a device under test (DuT) for our laser-based characterizations. The DuT is a hardware module implementing the Advanced Encryption Standard (AES) algorithm used, for example, for encrypting the secure channel between a SECRICOM host device and the SECRICOM SDM.

The AES algorithm

The The countermeasure against faults attacks implemented in this chip consists first in detecting errors and then reacting in case of detection. The error detection is done by using spatial duplication: the AES is executed twice in parallel and at each round, the results of the two instances are compared. If an error is detected, the reaction consists in blurring the erroneous cipher text with the scrambled value of the detected error. The error detection mechanism is described in the following section.

The figure below presents an overview of the ASIC AES architecture. We can see the two AES executed in parallel with the error detection system. The AES is interfaced with an AMBA APB bus to communicate with the external world. We used an FPGA for the interface between the APB bus and the PC used to control the AES and perform the different tests.

Fig. 6: Overview of the ASIC AES architecture

The error detection mechanism

As mentioned in the previous section, when an error is detected (XNOR operation between the two paths), the error detection mechanism scrambles the error value and then blurs the cipher text with the scrambled error. An error matrix is used whereby the error is spread across the rows and the columns as shown in the Figure 7. After that, the error matrix is XORed with the SubByte's results of the two data paths.

The Cross-ShiftRows operation

In addition to the above error detection mechanism, the ShiftRows operation is crossed between the two data paths. Half of the bits of each byte come from the other path and vice versa. This "crossing" operation is an additional security barrier. Indeed, if a fault is injected onto one of the paths, the fault is detected and propagated onto the two paths in parallel. With the Cross-ShiftRows, half of information is lost due to the transfer to the other path. Figure 8 illustrates the injection of a fault on the last round of the AES and its propagation throughout the two data paths.

Fig. 8: Propagation of a fault on the last AES round (with the countermeasures)

Fault injection using laser

The Laser test bench

To perform the different tests on the AES chip, we used the laser facility of the MicroPackS™ platform (MICROPACKS 2012). This laser is a YAG (Yttrium Aluminium Garnet) laser with three different wavelengths: green, infra-red and ultra-violet. We used the green wavelength (~532nm) with a 20x Mitutoyo lens. We obtained a spot size between 1µm and 150µm. With the larger spot size, we have an energy of approximately 15pJ per laser shoot.

The AES is interfaced with a control PC. When an encryption is launched, a trigger signal is sent to an FPGA synchronization board, which sends a shoot signal to the laser after a delay defined by the control PC.

This delay allows triggering the laser at different times during the encryption.

We put ourselves into two configurations when doing those tests. In the first configuration we adopted a "black box" approach where we ignored any of the implementation information we had on the DuT and tried to perform fault injections on the data path of the AES, in the 'classical one', with the objective of collecting the erroneous ciphertexts and doing differential cryptanalysis like in Giraud's or Piret's methods. In the second configuration, we used our knowledge of the implemented countermeasures, in a so-called "white box" approach, to try to circumvent the security mechanism by, say, trying to fault the detection mechanism itself.

Fault injection on the data path

The easiest way to generate errors and trigger the detection mechanism is to inject faults into the register of the SubByte module. Despite the propagation of the error onto the two data paths and the loss of half of the information due to the Cross-ShiftRows, we can always try to use the faulty ciphertexts in Giraud's DFA [START_REF]DFA on AES[END_REF]. The sine-qua-none condition for this attack is to generate mono-bit faults (i.e. errors on only one bit of the State matrix). The error matrix can be found with a simple XOR operation between the correct ciphertext and the faulty one. Figure 8 illustrates this kind of fault injection.

Fault injection on the detection mechanism

Another way of generating errors is to use the error detection mechanism itself. The DFA described in [START_REF] Piret | A Differential Fault Attack Technique Against SPN Structures, with Application to the AES and KHAZAD[END_REF] needs a fault injection before the last MixColumns operation (Round 9). With our countermeasures and a fault injection on the data path, the faulty ciphertexts are not exploitable for this attack. When looking closer at the Cross-ShiftRows operation, it appears that if the same fault is injected on the two data paths, the effect of the Cross-ShiftRows is "neutralized". Due to the dispersion of the lay-out of the two paths in all the chip's surface, it's very hard to inject the same fault into the two paths with a laser which has a local effect. The solution is to inject the fault directly into the error matrix. By doing so, we could propagate the same error on the two paths and the Cross-ShiftRows is neutralized and we do not have any loss of information. In the last round, as the errors are the same on the two data paths, no detection appears and we have a faulty ciphertext that could be used for DFA. Figure 9 depicts the propagation onto the two paths of an error injected directly into the error matrix at the 9 th round. Since the registers of the SubByte are dispersed across the ASIC, the first step of the characterization work was to physically localize those registers. To do so, we made a cartography using the laser. The surface of the chip was partitioned into 36 zones of size 150µm x 150µm each as shown in Figure 10. For each zone, 50 different encryptions were performed among which we looked for specific faulty ciphertexts.

Indeed, as shown in the Figure 8, if a fault is injected on one byte at the beginning of the last round, in the end, we have six times the same error, due to the propagation mechanism, and one different error at the position of the injected fault. In Figure 10, the colored the regions highlight those where specific faulty ciphertexts were obtained and which appeared to be more sensitive to this type of fault on the last round of encryption.

Results

Once the SubByte's registers have been localized, we first injected faults on the data path. In the "black box" approach, we tried to perform DFA as described in (GIRAUD 2005) using only the faulty and correct ciphertexts but in vain: the detection and error spreading mechanism proved to be efficient against such attacks.

However, when we used the knowledge of the implementation of the countermeasure (i.e. in a "white box" approach), especially the structure of the Cross-ShiftRows, we managed to recover a few bytes of the secret key of the AES. The complete knowledge of the Cross-ShiftRows is necessary because half of the information is lost in this operation and we need, for the attack, to keep all the information.

We also tried to inject faults directly on the error matrix in order to try another type of DFA [START_REF] Piret | A Differential Fault Attack Technique Against SPN Structures, with Application to the AES and KHAZAD[END_REF] where in theory two faulty ciphertexts are needed to recover 4 bytes of the secret key. To find all the 16 bytes of the key, we need to inject an error on one of the bytes of each column of the error matrix. Despite our effort, we couldn't inject a fault on the error matrix with the laser test bench. One of the reasons is that the error matrix is not implemented with registers but with logic gates. Thus it's very hard to synchronize the laser shoot with the ASIC's encryption precisely enough to target separately each column of the error matrix.

Discussion & Conclusion

In this paper, we have described how countermeasures implemented in a hardware implementation of the AES have been tested using a laser as fault injection means. We have seen that, in a "black box" approach, classical DFA techniques are inefficient against such countermeasures. However our characterization work has also shown that in a "white box" scenario, some bytes of the secret keys could be recovered. This has led us to the conclusion that the error propagation should have been truly random and independent from the generated errors (requiring the implementation of a True Random Number Generator on the chip). Moreover, another investigated attack path, trying to inject a fault on the error matrix itself, has proven to be unsuccessful illustrating the limits of current equipment with respect to current technologies.

Such characterization works have provided valuable design rules for implementing secure encryption AES modules like those used in the SDM of the SECRICOM project.

Fig. 1 :

 1 Fig. 1: Laser beam effect onto MOS structure Fig. 2: Transient current from charge collection after laser shoot

Figure 3

 3 Figure 3(a) shows an inverter with a low state on its input and the sensitive zone associated. The second inverter (b) on the figure represents the other case of localisation of a sensitive zone.

Fig. 3 :

 3 Fig. 3: CMOS inverter with a low state (a) or high state (b) on its input

  AES algorithm is a symmetric key cryptography standard established by the NIST(NIST 2001). This algorithm is based on a few transformations (SubBytes, ShiftRows, MixColumns, AddRoundKey) used iteratively in rounds (Figure5). In this paper, we focus on the 128-bit key version. This version processes data blocks of 128 bits, considered as matrices of 4x4 bytes called States, in ten rounds. The round keys (K1 to K10) used during each round, are calculated by a key expansion routine (not detailed in this paper). We denote M1 to M10 the AES States at the end of each round.

Fig. 5 :

 5 Fig. 5: The AES algorithm and its hardware implementation

Fig. 7 :

 7 Fig. 7: The Error Detection and Spreading Mechanism

Figure 9 :

 9 Figure 9: Injection of an error into the error matrix and its propagation

Fig. 10 :

 10 Fig. 10: Sensitive areas of the ASIC

Acknowledgements

This work was funded by SECRICOM project (EC FP7-SEC-2007 grant 218123). The research work of Cyril

Roscian is partly funded by the "Conseil Regional Provence-Alpes Côtes d'Azur".


Christophe Rizet

Cécilia Cruz

Matthieu De Lapparent

Martine Vromant

Cécilia Cruz

Cecilia Cruz

Carbon intensity of French shippers

Keywords: Green-House Gas, freight transport, shipper emissions. Résumé Gaz à effet de serre, transport de marchandises, établissement chargeur

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

Introduction : context and objective

Over one quarter of greenhouse gas (GHG) emissions in France comes from the transport sector and this share is growing. Furthermore, within the transport sector, the share of emissions from freight activity is growing and there is no sign of saturation (EMEP/EEA, 2010). Therefore deep changes are necessary to reach the target of dividing per 4 the anthropological GHG emissions, either through carbon efficiency improvements or by means of reduction of transport demand. The link between transport emission and climate change has been analysed in [START_REF] Schwanen | Scientific research about climate change mitigation in transport: a critical review[END_REF]; we are here focusing on freight carbon intensity. Carbon efficiency of freight transport and also the determinants of freight transport demand are not well understood. We claim that an in-depth knowledge of the factors that influence firms in their logistical choices as well as the factors that influence carbon efficiency would make freight GHG mitigation much easier. This lack of knowledge is due to several factors, including the theoretical complexity of the problem, the insufficiency of resources that have been made available for freight compared with passenger transport, and the inadequacy of the existing data. The French shipper survey developed in 2004 is an attempt to fill this gap in data. The survey methodology is discussed in [START_REF] Rizet | Energy consumption estimation with a shipper & transport chain survey[END_REF]; [START_REF] Guilbault | Transport and Logistics Demand: new input from large shippers surveys in France[END_REF] give an overview of the results.

Our objective in this paper is, using this survey, to estimate the level of CO 2 emissions and carbon efficiency of freight transport per shipper. After presenting the survey (section 2) and developing the way we estimated transport energy consumption and CO 2 emission per shipment (section 3), we present the two ways we used to estimate the yearly quantity of freight transport CO 2 emissions per shipper. First, we set up a typology of shippers, got 10 types of 'homogenous' shippers and estimated in each type an average level of emissions per shipper (section 4). In a second approach, we set up a model explaining the CO 2 emission per tonne shipped, using only the shipper characteristics as explaining variables. Multiplying this average CO 2 emission per tonne of the shipper per the yearly tonnage gives us the amount of emissions of this shipper for his shipments during that year (section 5). Then we compare the shipper emissions estimated through these two approaches and consider the variability of carbon efficiency among shippers (section 6).

Data : the 2004 French shipper survey

The two major characteristics of the French shipper survey are the description of the shippers' organizational features influencing the transport choices and the tracking of a selection of shipments from their departure from the plant up to their arrival to the consignee. In this survey, information is collected at three levels:

• Shipper establishment: after a few questions about the volume and structure of the company's ingoing and outgoing transport flows and its own fleet of vehicles, a face-to-face interview is administered to the logistics manager of the company. The questions concern the economic characteristics of the firm: production, distribution, storage practices, relationships with its customers and suppliers, and the management and communications systems it uses. A "transport" section is added to this description of the firm's industrial and logistical organization. It deals with the firm's relationships with carriers, terms of access to the various types of infrastructure, and how responsibility for transport is shared between the firm and its partners.

• Shipment level: at the end of the establishment questionnaire, the last 20 shipments are listed, of which 3 are randomly selected and then traced until they reach their final consignee. The consignment questionnaires are filled in either with the manager mentioned above or with the manager in charge of dispatching; they deal with the economic relationship between the shipper and the consignee and the terms of business between the two. This relationship between the shipper and his clients is described not only through the physical and economic characteristics of the shipment but also with information about the quantity of goods sent to this client per year (yearly tonnage to this consignee), the yearly number of shipments to this consignee, the split of responsibilities with regard to transport organization and the contractual allocation of transport costs and associated services. The first information required in order to reconstruct transport chains is also collected at this level, with the identification of the consignee.

• Transport operator and journey link: for this part of the survey, questionnaires are administered by telephone.

Questions deal with the economic characteristics of the operator, including the operator to whom it handed the shipment and with the characteristics of the transport leg: mode and vehicle type, load, etc. The next operator is in turn questioned up to the final consignee and the description of the whole transport chain. The transport chain is split into as many legs as there are changes to another vehicle.

The transport chain is therefore reconstructed throughout Western Europe and includes an interview with the consignee. For a shipment which travels beyond the limit of Western Europe, the transport chain is surveyed only until the first transfer point after the frontier of Western Europe. Furthermore, the 2004 shipper survey has been adapted to enable the quantification of energy consumed and CO 2 emitted, and to relate them to the determinants of freight transport demand [START_REF] Rizet | Energy consumption estimation with a shipper & transport chain survey[END_REF]. Primary and final energy consumption have been computed per leg, as the product of distance per vehicle energy consumption (for trucks: litres per km), divided per the load and multiplied per the shipment weight; empty backhauling is estimated using the results of the French national freight vehicles survey. Energy consumption is converted in grams of oil equivalent (goe) and shipment energy consumption is the sum of energy per leg. Both Tank to Wheel (TW) and Well to Wheel (WW) CO 2 emissions are then computed as the product of energy per an emission factor. Dividing CO 2 per the tonnekilometre of the shipment gives the carbon efficiency of the shipment, in kg CO 2 /tkm. The method to compute CO 2 per shipment is described in [START_REF] Rizet | Impact of logistical choices on freight transport carbon efficiency WCTR 2013 Rio[END_REF]. In this paper we present the ways we estimated the yearly quantity of freight transport CO 2 emissions per shipper. This was done with two approaches: First we set up a typology of shippers, in 10 types of 'homogenous' shippers and estimated the average level of CO 2 emissions in each type; in the second approach, we modelized the carbon intensity per tonne shipped and multiplied this per the yearly tonnage shipped by a shipper.

Estimating shipper yearly emissions from a typology of shippers

In the survey, the sample of shipments is obtained by a first sampling among shippers (about 3000 surveyed shippers) and a second sampling of shipments among the surveyed shippers. Only three shipments are surveyed per shipper on average, which does not allow estimating directly the volume of emissions of this shipper. So we grouped a significant number of shippers having 'profiles' as close as possible, to obtain a number of shipments which allows to estimate the volume of emissions of this type of shippers.

A typology of shipper establishments

Our objective thus is to establish a typology of establishments using the variables which characterize the volume of emissions of these shippers. After different tries, we finally used a typology realized by 'analysis of the multiple correspondences' (AMC) and a 'non hierarchical classification', in 10 types [START_REF] Clausen | Applied correspondence analysis[END_REF]. We used the following variables:

• The number of employees • The annual number of shipments • The distribution of all the shipments according to their weight, • The distribution of all the shipments according to their destination, • The distribution of all tons shipped according to the main transport mode.

AMC allows realizing a typology based on the resemblance. Two establishments will be close all the more as they possess a large number of modalities (variables) together. The observation also concerns the association between modalities. We did it with the SAS CORRESP procedure. The output of this procedure is a set of new variables which correspond to the coordinates of the individuals (establishments) and the variables in a multidimensional space. This space is defined by more or less numerous axes according to the number of variables or the quality of the analysis. The quality of the explanation of the cloud of points will be a function on behalf of slowness of this cloud which is explained by axes. For the classification we used 6 axes, what brings the part of slowness explained to 51 %. The typology is then realized by the method of classification of dynamic clouds applied to the coordinates of the establishments on each of these axes. The next step is to group shippers by types, according to their positioning on the axes. Two establishments close according to these axes will be grouped in the same type. The dynamic clouds method is then applied to the coordinates of shippers on each of these axes. Shippers are grouped in a limited number of types in such a way that the establishments of the same group are alike from the point of view of variables used for the analysis. We chose to limit to 10 the number of types.

Starting from 10 cores and by successive iterations, the procedure supplies a set of types by grouping the closest establishments there. The SAS CLUSTER and FASTCLUS procedures were used to define the types from the coordinates stemming from the AMC. We so obtained 10 types of shippers with more than 200 shippers in each type, plus a group of 120 shippers with missing variables, classified as 'without type'.

Estimating the quantity of freight transport CO 2 emissions per shipper of each type

For each type of establishments so obtained we calculate then (cf. table 1 below):

• The total quantity of CO 2 emissions of freight transport as the weighted sum of the emissions of the shipments of this type;

• The average emissions of an establishment of this type, by dividing the total emissions of the type by the weighted population of establishments.

Levels of emissions are very variable among types. Biggest emitters are in Type 2 with on average 4764 tonnes of CO 2 per year for the transportation of their shipments; the lowest emitters are in type 1, with an average of 79 tonnes per establishment and per year. To appreciate the dispersal of these values of average emissions per shipper according to the type, we calculated for each type, the ratio between the average emissions per shipper of this type to the emissions of the least emitting type shipper (i.e. type 1): according to this ratio (last column in table 1), an establishment of type 2 emits on average 56 times more that an establishment of type 1. The various types do not all differ significantly by their average level of emissions: types 5 and 9 on the one hand and types 7 and 10 on the other hand seem to have close levels of yearly emissions.

Precision of the estimated quantity of CO 2 per shipper

To calculate the precision of this estimator, we used the 'bootstrap' method: in each type we randomly sampled important groups of establishments and, on each of these sub-samples we calculated a new value of the CO 2 by establishment with the previous method. The distribution of these indicators from the various sub-samples of the same type allows us to calculate the precision of the estimation of the CO 2 per shipper of this type. The SAS SURVEYSELECT procedure allows making this poll of sub-samples and so facilitates this method of 'bootstrap'. Figure 1 shows the distribution of the average CO 2 by establishment so estimated in the various subsamples of type 1. The average values and the standard deviation of these estimators resulting from the bootstrap, enable to calculate the limits of the 95 % confidence interval. These limits are shown for each type in figure 2 here under. Figure 2 shows that only the type 2, with an average of 4770 t of CO 2 per shipper and per year clearly differs from other types for the volume of emissions. Types 1 and 8, who have a low level of emissions, are very close one of the other, with an average of respectively 80 and 87 t. and confidence intervals which widely recover. Confidence intervals of the other types also recover: types 3, 4, 5, and 9 on the one hand and types 6, 7 and 10 on the other hand.

Estimating shipper yearly emissions with a model of carbon intensity per tonne shipped

We want to develop and calibrate a model linking the CO 2 emitted by the transport to the shipper characteristics. If we disregard the estimation of the levels of emissions by type of shippers, transport CO 2 emissions are known for shipments. So, various models were tested to explain the carbon intensity of a shipment by the characteristic of the shipper. As far as the only explanatory variables are the shipper characteristics, such a model enables to estimate the average carbon intensity of a shipment from this shipper and, by multiplying this carbon intensity per the number of shipments, to estimate the total emissions of this shipper per year. In practice, we modelled the carbon intensity per tonne shipped rather than per shipment.

A model of carbon intensity per tonne shipped

Two quantities are known for all shippers and could allow estimating their freight transport emission: the yearly number of shipments and yearly tonnage shipped. After various tries, we decided to model the carbon intensity per tonne, rather that the intensity per shipment, because the estimate is more accurate. Various forms of models were tested, in particular combinations of logarithms and levels.

The potential explanatory variables are the characteristics of the shipper from the survey, in order to be able to estimate an average level of carbon intensity per tonne for each shipper, independently of the characteristics of the shipment. In these 300 available variables, about 15 were found explanatory without counting too many missing It is based on 7 573 observed shipments with all the explanatory variables; the coefficient of correlation is 0.40 and 14 variables are explanatory at the level of 5 %. These 15 variables are the characteristics: -Of the distribution of the shipper shipments according to their weight (a proxy for the shipment weight); -Of the distribution of the shipper shipments according to their geographical destination (a proxy for distance); -Of how much this shipper is using the various transport modes (a proxy for the mode) and -Of the shipper 'business sector' (a proxy for the type of goods shipped). Estimated carbon intensity per tonne varies from about ten kilograms to some 3 tonnes with a great majority below one tonne.

Precision of the estimated carbon intensity per tonne

The 95% confidence interval of carbon intensity was calculated with SAS, using the option ' lclm ' = lower limit and ' uclm ' = upper limit. In figure 3 below, every establishment is represented by three points which have for abscissa the estimated CO 2 /tonne of the shipper and, for ordered, the same value of carbon intensity per tonne shipped by this shipper, as well as the lower and higher limits of the 95% confidence interval of this estimated value. To allow a better legibility, the figure 3 The limits of the confidence interval are about 20 % above and below the estimated value. This confidence interval varies from an establishment to the other one, but this variability is limited itself: there is thus an important difference between establishments on the left of the figure, which emit less than 100 kg by transported tonne and the establishments on the right, which emit more than 900 kg by tonne.

The limits of the confidence interval are about 20 % above and below the estimated value. This confidence interval varies from an establishment to the other one, but this variability is limited itself: there is thus an important difference between establishments on the left of the figure, which emit less than 100 kg by transported tonne and the establishments on the right, which emit more than 900 kg by tonne.

Carbon intensity per shipper Carbon intensity per shipper

The CO 2 emitted by a shipper per year is calculated as the product of the CO 2 per tonne, estimated by shipper according to the method described above, by the yearly tonnage shipped by this establishment. The yearly tonnage of the shipper is known from the survey so the limits and the confidence interval for the CO 2 per shipper are calculated as the product of the limits of the confidence interval of the CO 2 by tonne multiplied by the yearly tonnage of the shipper.

The CO 2 emitted by a shipper per year is calculated as the product of the CO 2 per tonne, estimated by shipper according to the method described above, by the yearly tonnage shipped by this establishment. The yearly tonnage of the shipper is known from the survey so the limits and the confidence interval for the CO 2 per shipper are calculated as the product of the limits of the confidence interval of the CO 2 by tonne multiplied by the yearly tonnage of the shipper. According to this approach, the CO 2 emitted by a shipper in one year for the transport of its shipments varies from approximately one tonne (on the left of figure 4) up to more than 10 000 tonnes (on the right). The confidence intervals are important in absolute value but, considering the dispersal of these estimated values of the CO 2 per shipper, it is possible to define groups of establishments clearly differentiated from the point of view of emissions, with little overlapping of the confidence intervals.

Precision of the estimated carbon intensity per shipper Precision of the estimated carbon intensity per shipper

According to this approach, the CO 2 emitted by a shipper in one year for the transport of its shipments varies from approximately one tonne (on the left of figure 4) up to more than 10 000 tonnes (on the right). The confidence intervals are important in absolute value but, considering the dispersal of these estimated values of the CO 2 per shipper, it is possible to define groups of establishments clearly differentiated from the point of view of emissions, with little overlapping of the confidence intervals.

Comparison of carbon intensity per shipper estimated by the two methods Comparison of carbon intensity per shipper estimated by the two methods

We have now two values of the CO 2 emitted by shipper for the transport of its yearly shipments, one established by the typology and the other one by the model. Figure 5 below compares these indicators by type of shippers, for those shippers for which a level of CO 2 has been estimated by the two ways. The biggest emitters (shippers of type 2), with 4700 t per year emit 60 times more that those of type 1 (78 t). A CO 2 mitigation policy aiming at

We have now two values of the CO 2 emitted by shipper for the transport of its yearly shipments, one established by the typology and the other one by the model. Figure 5 below compares these indicators by type of shippers, for those shippers for which a level of CO 2 has been estimated by the two ways. The biggest emitters (shippers of type 2), with 4700 t per year emit 60 times more that those of type 1 (78 t). A CO 2 mitigation policy aiming at shippers should clearly focus on these 'big emitters'. On average the CO 2 emitted by a shipper is a little more than 620 t a year.

Carbon intensity, carbon efficiency and the impact of a carbon tax on a shipper competitiveness

To quantify CO 2 emitted by shippers is important for assessing the potential impact of a carbon tax policy. If we consider that such a carbon tax would be distributed on all the tonnes sold by this company, then the main indicator for the impact on shipper competitiveness should be the shipper CO 2 per tonne, i.e. its carbon efficiency. If we now consider that the amount of carbon taxes will be used to reduce the cost of labour (health, unemployment, and pension), in order to improve the competitiveness of the national workers, then the benefit of such a policy for a company would be proportional to its number of employees: the more employees, the more cost reduction.

Our 10 types of shippers not only have different carbon intensity per shipper and per employee but they also correspond to different carbon efficiencies per tonne (from 12 up to 208 kg CO 2 /tonne) and per tkm. Carbon efficiency per tkm is less variable (from 0.046 up to 0.678 kg CO 2 /tkm). By employee, the average level of emission, is 16.5 t CO 2 per year and varies from 4.3 t for type 10 up to 27.2 t for type 9 that is a 1 to 6 ratio. Indeed, these four indicators give rather different views of the burden of a freight transport carbon tax among the different types of shippers.

• A transport carbon tax might be supported by the shippers according to their level of emissions. With a tax of 100 € per tCO 2 for example, the tax per shipper would range from 7 900 € for shippers of type 1 up to 476 400 € for shippers of type 2. • These shippers have very different levels of sales or shipments; related to their tonnage shipped a tax of 100 € per tCO 2 would then range from 1.18 €/sold tonne for shippers of type 2 up to 20.8 €/sold tonne for shippers of type 10.

• If the product of this carbon tax is used to reduce the cost of labour, the benefiters would then be the big employers. The number of employees is then an indicator of the cost reduction for the shipper of such a policy.

• The ratio tCO 2 emitted per employee is then a good indicator of the real cost for a given shipper or type of shipper. Multiplied per the ratio [tax (in €/tCO 2 ) / labour cost reduction (in €/employee), it gives an idea of the impact of such a policy for this type of shippers. If the product of these 2 ratios equal one, the benefit is as important as the cost. If this product is higher than one then the cost is higher than the benefit; and if this product is smaller than one, then the reduction of cost is more important than the tax paid. We don't know the ratio 'tax per tonne/reduction of cost per employee' but we can make the hypothesis that this ratio will be identical for all shippers. Then the comparison of the tCO 2 /employee ratio is an indicator of the split of the burden among the different shippers. The winners of such a policy would then be the shippers of type 1, 3 and 8; the shippers of type 10 would be the main loosers.

• Carbon efficiency, gives another view of the comparison of emissions between very heterogeneous shippers.

Type 10 for example has e very low carbon efficiency per tonne and per tkm.

Conclusion :

CO 2 emissions per shipper has been estimated with two different approaches, using the French shipper survey data: First we set up a typology of the shippers, in 10 types of 'homogenous' establishments and estimated in each type an average level of emissions per shipper. In a second approach, carbon intensity per tonne shipped has been modeled, using only the shipper characteristics as explaining variables. Then for each shipper, multiplying this average CO 2 emission per tonne per its yearly tonnage shipped gave the amount of emissions for its shipments during that year. The two approaches gave reasonably similar results. This in depth knowledge of the shipper CO 2 emissions can be used to compare the carbon efficiency of different types of shippers and the impact of a carbon tax applied to freight transport on shipper's competitiveness. Four indicators of carbon performance have been quantified that give very different results for our ten types of shippers:

Carbon intensity per shipper and per year varies from 78 tCO 2 , for the 'lower emitter type' up to 4700 t for the 'upper' type. The average is 620 tCO 2 per shipper and per year and a shipper of the 'upper' type emits 60 times more CO 2 than the 'lower' ones. The ten types of shippers also have very different carbon intensities per tonne: from 12 up to 208 kgCO 2 /tonne (a 1 to 17 ratio) with an average of 46.7 kgCO 2 /tonne. Per tkm, carbon intensity varies according to the shipper types, from 46 up to 678 gCO 2 /tkm, i.e. a 1 to 15 ratio, with an average carbon intensity of 108 gCO 2 /tkm. Finally, carbon intensity per employee could be important, if the carbon taxes receipt would be used to reduce the indirect cost of labour: shippers would pay taxes proportionally to their emissions and would benefit reduction of labour costs, according to their number of employees. Then group of shippers' n° 1, with a ratio of 4.3 tCO 2 per employee per year would be the winners while shippers of group 10, which have a ratio of 27.2 t/employee, would be the main losers of such a policy

Figure 1 :

 1 Figure 1 : Distribution of the estimators of CO 2 per shipper of type 1, in tCO 2 / shipper / year

Figure 2 :

 2 Figure 2 : Precision of the CO 2 estimator per shipper according to his type in tCO 2 /shipper (the blue bars correspond to the minimum of the confidence interval and the red ones to the maximum)

  . The class variables have been transformed in dummies. The final model is: Log (CO 2 /tonne) = a1* x1 + a2* x2 + a3* x3 + ...+ a15* x15 Where x1, x2, x3... are the shipper characteristics.

Figure 3 :

 3 Figure 3: Carbon intensity per tonne shipped and the 95% confidence intervals Figure 3: Carbon intensity per tonne shipped and the 95% confidence intervals

Figure 4 :

 4 Figure 4 : Estimated tonnes of CO 2 per shipper and the 95% confidence intervals Figure 4 : Estimated tonnes of CO 2 per shipper and the 95% confidence intervals

Figure 5 :

 5 Figure 5 : A comparison of shipper carbon intensity estimated by the two ways

Table 1 :

 1 The 10 types of shippers and their CO 2 Emissions

	Type number	Sample size (number of observed shippers)	Population (weighted number of shippers)	Total emissions of the type (thousands tCO 2 WW)	Average emissions / shipper (tCO 2 WW)	Ratio
	0 (no type * )	120	1 869	1 970	1054	13,3
	1	287	13 731	1 089	79	1
	2	265	2 228	10 614	4 764	60.3
	3	227	4 999	1 261	252	3.2
	4	399	7 260	4 632	638	8.1
	5	300	4 626	1 860	402	5.1
	6	324	3 340	5 189	1 554	19.7
	7	246	3 810	4 612	1 210	15.3
	8	236	15 127	1306	86	1.1
	9	261	4 897	2 421	494	6.2
	10	214	7 335	8 031	1 095	13.8
	TOTAL	2879	69 222	42 985	621

*

: Shippers with missing data

  focuses on carbon intensities lower than one tonne.

	kg CO 2 /t.
	kgCO 2 per shipped tonne

Table 2 :

 2 carbon intensity and carbon efficiency per type of shipper

Acknowledgements

The results presented here are part of an ongoing Ifsttar research program funded by the French Ministry of transport (DGITM); we have used the SAS software for data analysis and figures.


Liu Bo

Franck David

Christophe Coquelet

Renaud Gicquel

Philippe Rivière

Performances and cost assessment of an air cooled two-stage rankine cycle for large power plants operating with different working fluids

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


R Tayade

A Pasko

C Serpico

F Mazaleyrat

M Lobue

Magnetization reversal in exchange-spring bilayer system under circularly polarized microwave field

Keywords: Micromagnetics, Magnetic switching, Microwave-assisted reversal

Microwave assisted magnetization reversal are studied in the bulk bilayer exchange coupled system. We investigate the nonlinear magnetization reversal dynamics in a perpendicular exchange spring media using Landau-Lifshitz equation. In the limit of the infinite thickness of the system, the propagation field leads the reversal of the system. The reduction of the switching field and the magnetization profile in the extended system are studied numerically. The possibility to study the dynamics analytically is discussed and an approximation where two P-modes are coupled by an interaction field is presented. The ansatz used for the interaction field is validated by comparison with the numerical results. This approach is shown to be equivalent to two exchange coupled macrospins.

I. INTRODUCTION

ith the tremendous demand in the electronic storage technology, new strategies are being continuously tested to overcome the so called "recording trilemma" [START_REF] Weller | Thermal effect limits in ultrahigh-density magnetic recording[END_REF]. This fastpaced developing technology requires understanding of the fundamental dynamical magnetization reversal processes [2] [3]. One of the current challenges is to reduce the writing field in the recording devices as the use of high perpendicularanisotropy materials requires very intense applied fields for switching the media grains. One of the projected solutions is to use the heat-assisted magnetic recording (HAMR) or thermally assisted recording (TAR) technique [START_REF] Moser | Magnetic recording: advancing into the future[END_REF]. Thirion et al. [START_REF] Thirion | Switching of magentization by nonlinear resonance studied in single nanoparticles[END_REF] proposed a method to overcome high fields required to reverse the spins in the high density media without compromising the thermal stability of the devices by introducing a radio-frequency (RF) field along with a constant applied field [START_REF] Tanaka | Micromagnetic study of Microwave-Assisted Magnetization Reversals of Exchange-Coupled Composite Nanopillars[END_REF]. Suess et al. [START_REF] Suess | Exchangecoupled perpendicular media[END_REF] reviewed the use of exchange spring media for magnetic recording and moreover they showed that it is possible for perpendicular recording media to achieve the required thermal stability at very high densities.

In this work we focus on the application of microwave assisted switching to perpendicular exchange spring media (PESM). Many works have been devoted to show the reduction in the switching field in PESM compared to single phase high anisotropic media [START_REF] Bashir | Microwave-Assisted Magnetization Reversal in Exchange Spring Media[END_REF], [START_REF] Igarashi | Mechanism of microwave assisted magnetic switching[END_REF], [START_REF] Boone | Microwave-Assisted Magnetic Reversal in Perpendicular Media[END_REF]. However, there is very limited literature to explain analytically the mechanisms of switching and reversal in this class of materials [START_REF] Suess | Micromagnetics of exchange spring media: Optimization and limits[END_REF]. The analytical description of magnetic behavior in PESM was studied by Asti et al. [START_REF] Asti | Magnetic phase diagram and demagnetization processes in perpendicular exchange-spring multilayers[END_REF] using the static micromagnetic approach. The magnetization dynamics was analytically studied by Bertotti et al. [START_REF] Bertotti | Nonlinear-dynamical-system approach to microwave-assisted magnetization dynamics (invited)[END_REF], [START_REF] Li | Analysis of Microwave Assisted Magnetization Switching in Magnetic Material[END_REF] for uniaxial anisotropic media. This work was limited to uniformly magnetized systems.

The challenge represented by the non linear nature of Landau Lifshitz (LL) equation, and by the non-uniformity of the problem, due to the interface between the hard and soft magnetic phases, makes a pure analytical solution hardly achievable. However, here we show that the numerical solution of the problem makes possible, under rather general condition, to describe the dynamics using the same analytical results obtained for the decoupled layers with an additional interaction field mimicking the effect of the interface induced nonuniformity.

II. MODEL

We consider a system composed of two layers lying on the xy plane, one magnetically hard and the other soft. Both layers exhibit perpendicular uniaxial anisotropy directed along zaxis. The system is infinite in the xy plane and due to this geometry the magnetization is uniform along the plane. Thus the problem reduces to work out the distribution of magnetization vector M= M(z) with respect to z, is the local magnetization vector. Calling z 0 the coordinate of the interface between the soft and the hard layer we have that |M(z)| = M s (z). The in-plane component of the magnetization is given by M  = M 2

x +M 2 y . The magnetocrystalline anisotropy is K(z) = K 1 for z < z 0 , and K(z) = K 2 otherwise. The shape and magnetocrystalline anisotropies can be merged in a single effective anisotropy constant, K u (z) = K(z) - 0 M s 2 (z)/2. The two phases are assumed to be exchange coupled at the interface.

The magnetization dynamics are given by LL equation,

) ( ) ( ) 1 ( 2 eff s eff z M t                   . ( 1
)
where γ is the gyromagnetic ratio associated with electron spin, α is the damping constant, and H eff is the effective field, , ) (

) ( 2 ) ( 2 2 2 0 z l z M z K t ex s u a eff         z e  . ( 2
)
where

l ex = (2A/ 0 M 2 s (z)) 1/2 is the characteristic exchange length.
H a (t) is the total applied field which is the sum of a DC magnetic field, H az applied perpendicular to the plane of the media and a circularly polarized microwave field of magnitude H a rotating with angular frequency  = 2πf, W applied in the plane.

). sin (cos ) (

y x a z az a t t H H t e e e        (3)
It is known [START_REF] Asti | Magnetic phase diagram and demagnetization processes in perpendicular exchange-spring multilayers[END_REF] that there are two fields controlling the process of magnetic reversal in spring systems, the nucleation field H N and propagation field H P . For exchange coupled multilayer systems it has been shown that the leading field in the reversal process, for thin systems is H N and H P in the limit of infinite thickness of the layers. In this case, the analytical expression for the propagation field in static conditions is [START_REF] Asti | Magnetic phase diagram and demagnetization processes in perpendicular exchange-spring multilayers[END_REF],

.

1 1 2 1 2 0 1          u u s u P K K M K H  (4)
To study numerically the dynamical reversal of this system we have developed a Fortran code solving LL equation using the finite difference technique [START_REF] Serpico | Numerical technique for integration of the Landau--Lifshitz equation[END_REF]. Here we focus our attention on a thick system where the propagation field leads the reversal process. Actually the infinite thickness limit is easily reached when the total width of the system is few tens of its domain wall.

The switching process is studied keeping H a and  constant and by slowly changing the DC field between H az = -1.6 × 10 6 A/m and H az = 1.6 × 10 6 A/m. The effect of eddy currents on the uniformity of the applied field is neglected.

We solve the equation of motion for the bilayer system using the micromagnetic parameters of a widely used PESM, FePt/Fe system. The parameters used are as follows: exchange constant A 1 =A 2 -11 J/m, anisotropy constants K 1 =2 × 10 6 J/m , 2 = 0.48 × 10 5 J/m 3 . We performed calculations assuming same saturation for both layers M s = 1.55 × 10 6 A/m (this overestimates M s for hard layer but allows a simpler analytical formulation in term of reduced units). Damping parameter is taken to be α = 0.01. Soft and hard layers have equal thicknesses of 100 nm. The magnitude of the in-plane microwave field is kept constant, H a = 8 × 10 4 A/m. Simulations are performed for several microwave frequencies ranging 1 GHz  f  20 GHz.

III. RESULTS AND DISCUSSION

A first and worth discussing result of the simulations is that, in a frequency interval up to more than 10 GHz, the average magnetization of the system, taken along z, <M z >, behaves like a P-mode. A P-mode is a stable solution of the dynamic equation for the system subjected to the microwave field, which results in uniform precession of the magnetization around anisotropy axis, in synchronization with microwave field. Due to non-uniform magnetization in the bilayer system, the stable solution is given by the average magnetization <M z >, which is constant in time forming an angle  with respect to z axis whereas < M  > is rotating with angular frequency , and with a phase lag  with respect to the rotating component of the applied field H a . We call this stationary state a global P-mode. In Fig. 1 we plot the normalized average z-component of the magnetization against the DC magnetic field showing the ascending branch of a switching loop with -1.6 × 10 6  H a  1.6 × 10 6 A/m. The Pmode behavior is apparent up to f = 10 GHz. Beyond this limit, the average magnetization shows oscillations whose behavior is much more similar to a quasi periodic mode (Qmode, using terminology from [START_REF] Bertotti | Nonlinear Magnetization Dynamics under Circularly Polarized Field[END_REF] [START_REF] Mayergoyz | Nonlinear magnetization dynamics in nanosystems[END_REF]. Another feature apparent in Fig. 1 is the switching field reduction as a function of the microwave frequency (see inset in Fig. 1). This microwave assisted switching is similar to the one reported in [START_REF] Bertotti | Nonlinear-dynamical-system approach to microwave-assisted magnetization dynamics (invited)[END_REF] for uniformly magnetized systems. It is worth noting that 6 GHz curve is switching at a very small positive DC field. This results opens a possibility to investigate zero field directional field microwave assisted switching.

A deeper analysis of the simulation results gives us the opportunity to go beyond the global P-mode concept. Fig. 2 shows the angle θ(z) between perpendicular component of magnetization M and H az as a function of z for different values of H az . The system is subjected to a constant in-plane circularly polarized field with f = 6 GHz. The figure shows two uniform regions separated by a thin transition region. The thickness of the latter stays constant (~30 nm) all along the magnetization curve and its contribution to the average magnetization can be assumed to be negligible. A plot of the lag angle (z) between M  (z) and H a gives very similar results. Thus we can say that the global P-mode is nothing but the average of two individual P-modes, corresponding respectively to the soft and the hard layer. Treating with P-mode solutions, it is possible to greatly simplify expressions by changing the system frame from the laboratory frame to a frame in which the system is rotating around z axis with the same angular frequency  of the microwave field. In this frame, the applied field becomes constant and is given by H a = H az e z + H ax e x and the LL equation ( 1) is written as, ),

2 eff s z eff z M t                               e ( ) ( ) 1 ( ) 1 ( 2
The magnetization M is fully defined by its angle θ with respect to z-axis and the angle  with respect to x (i.e. to H a ). From [START_REF] Bertotti | Nonlinear Magnetization Dynamics under Circularly Polarized Field[END_REF], and [START_REF] Mayergoyz | Nonlinear magnetization dynamics in nanosystems[END_REF] we know that a P-mode maps, for given , the values of M z and of the lag angle  (namely of M x , and M y ) onto the values of the DC field H z , and on the in-plane rotating field amplitude H  through the following expressions: From Eq. ( 6), [START_REF] Suess | Exchangecoupled perpendicular media[END_REF], and (8), we can easily verify that the field we are actually applying, H a would not produce the observed P-mode when applied to the uncoupled layers. From this we can suppose the observed P-mode to be generated by a field H i (subscript i = 1, 2 refer to the P-mode we observe in i layer) being the sum of the applied field and of an interaction field H j, exc representing the interface exchange coupling with the other layer, .

z s u z M M K H          0 0 2     (6)   2 2 2 2 2 2 2 0 z s s M M M H                 (7)
,exc

j a i H H H   (9)
We use a standard ansatz for the exchange field (see [START_REF] Suess | Micromagnetics of exchange spring media: Optimization and limits[END_REF]) with H j,exc = C M j , where C is a constant to be determined. Assuming that,

        L z z if t z z if t t z 2 2 1 1 ) ( 0 ) ( ) , (    (10)
where, L is the total thickness of the system and z 1  z  z 2 , is the region where the transition layer is located. We define |M 1,2 (t)| = M s1,2 . The interaction due to the layer is taken into account by an energy term - 0 CM 1 M 2 , so that the effective field in the two spatially uniform layers is given by

1 , 2 2 , 1 2 2 , 1 0 2 , 1 2 , 1 1 2    C e M M K z z s a eff              (11)
The P-mode in the two coupled layers is given by the equations 0

) ( ) 1 ( 2 , 1 2 , 1 2 , 1 2 , 1 2 2 , 1 2 , 1                eff s z eff M          e (12)
The four unknowns of the problem, specifically the three components of H i and C, are fully determined by equations ( 6), ( 7), ( 9) and by the constraint, H j,exc  = C M s . Fig. 3 Representation of exchange energy versus the angle between M 1 and M 2 defined in Eq. [START_REF] Boone | Microwave-Assisted Magnetic Reversal in Perpendicular Media[END_REF].

Using our numerical results we calculated, C = 0.1. On the other hand, we can also estimate C by calculating the exchange energy of the bilayer system. The exchange energy term is calculated by integrating (M) 2 along the sample thickness. This integration includes the transition region, which, of course, will give the main contribution to the total exchange energy of the system. Exchange energy is plotted against the angle between the magnetizations of the two uniform regions (i.e. the two P-modes), respectively M 1 and M 2 , Fig. 3. The slope of the curve gives the exchange coefficient C. Linear fitting of the curve is shown by the red line in the figure .  The points in the figure which do not follow the linear behavior correspond to the region around switching field. This approach is equivalent to describe the system as two exchange coupled macrospin with different uniaxial anisotropy. This model is numerically solved under the same parametric conditions used for the bilayer system with the exchange constant derived from the solutions of bilayer system. In Fig. 4 we show that <M z >(H az ) loops for the extended bilayer and of the macrospin system are in very good agreement.

Further investigations will be devoted to apply the bifurcation analysis presented in [START_REF] Bertotti | Nonlinear-dynamical-system approach to microwave-assisted magnetization dynamics (invited)[END_REF] to the case of two coupled P-modes. In this case, the phase portrait analysis discussed in [START_REF] Bertotti | Nonlinear-dynamical-system approach to microwave-assisted magnetization dynamics (invited)[END_REF] can be used to describe the dynamics of each single P-mode substituting the field of Eq. ( 9) to the standard applied field. Fig. 4 Average magnetization vs DC field. The loop in red (solid line) is the magnetization behavior for bilayer system obtained from the numerical simulations and the loop in blue (with filled circles) is the magnetization behavior from the macrospin magnetization model.

Fig. 1 .

 1 Fig.1. Spatial average of M z as a function of frequency. Only the branches of the loop when DC field is varied from -1.6 × 10 6 to 1.6 × 10 6 A/m, is shown here for frequencies of 1,6,10 and 15 GHz. Variation of switching field in response to the microwave frequency (line with triangles) is shown in the inset with the solid straight line corresponding to the switching field in absence of microwave field.

Fig. 2

 2 Fig.2 Spatial distribution of θ(z) for various DC field values at AC field frequency 6GHz. From z=0-100 represents the soft phase and z=100-200 represents hard phase. The black circles at θ=π represents the initial state of the system and the down triangle represents the final state of the system after complete switching.

ACKNOWLEDGMENT

The authors wish to thank Giorgio Bertotti for many helpful suggestions. This work was funded by the ANR project MICROMANIP and by the Fondation de Cooperation Scientifique Campus Paris Saclay as a part of the foreign guest project of DIGITEO N• 2012-XXD. This work was partially supported by MIUR-PRIN Project No.2010ECA8P3 "DyNanoMag."