Julie Mougin

Geraud Cubizolles

Anne Hauch

Jari Pennanen

Joaquin Alvarez

Sergii Pylypko

Marc Potron

Bastien Marquillier

Stephane Hody

Gianna Cesareo

Sabina Fiorot

Guillermo Perez

Development of an efficient rSOC based renewable energy storage system

While renewable energies locally installed in buildings are increasing; there is a potential or even a need to store them to promote their self-consumption. Being efficient, modular and scalable, reversible Solid Oxide Cell (rSOC) technology can play a key role. Works undertaken in the frame of the REFLEX European project are presented, starting with the optimization of each individual component, ie. cells, stacks and BoP components such as power electronics. With the support of modelling activities, the rSOC based system, made of 3 modules of 4 stacks each, coupled with a battery storage has been designed in such a way to ensure achieving maximum efficiency, while still operating the system in safe and lifetime-optimal conditions. The site where the system will be installed has been prepared. Furthermore, techno-economic simulations have been carried out to evaluate market driven requirements on system sizing and costs to be competitive with other storage solutions.

Introduction

According to the European Green Deal and Europe's clean energy transition, renewable electricity is expected to decarbonise a large share of the EU energy consumption by 2050 [START_REF]The European Green Deal[END_REF], but not all of it. Hydrogen has a strong potential to bridge some of this gap, as a vector for renewable energy storage, alongside batteries (2).

On 25 February 2015, the European Commission adopted "A Framework Strategy for a Resilient Energy Union with a Forward-Looking Climate Change Policy" (3). Resilient energy systems need to be developed and deployed, with a central aspect being the energy storage of rapidly increasing intermittent and low predictable renewable energies. According to IEA [START_REF]Energy Technology Perspectives 2016 -Towards Sustainable Urban Energy Systems, IEA[END_REF], as the share of the world's population living in cities rises, ambitious action in urban areas can be instrumental in achieving long term sustainability of the global energy system -including the carbon emission reductions required to meet the climate goals defined at COP21 in Paris. Indeed, buildings and eco-districts are increasingly becoming energy producers, with growing amounts of PV panels and wind power installed. Optimizing the self-consumption of this locally produced energy becomes a priority issue for the energy transition in order to avoid putting strain on the transmission and distribution network. Storing this energy to restore it locally, while having the option to produce locally with the same equipment combined heat and power (potentially also from methane or biomethane and not only H2), will ensure that the electricity supply matches demand at all times. The reversible Solid Oxide Cell (rSOC) technology allows this. Indeed, it is able to operate either in SOEC (solid oxide electrolysis) mode to store excess electricity to produce hydrogen or in SOFC (solid oxide fuel cell) mode in order to produce electricity (and heat) again, from hydrogen or any other fuel locally available when energy needs exceed local production. The rSOC system is advantageously completed with an electrochemical storage solution allowing fast response to the electrical energy needs.

In this context, the REFLEX European project will demonstrate, in-field, the high power-to-power (P2P) round-trip efficiency of this rSOC based technology (as compared to other H2 based solutions) and its flexibility and durability in dynamic operation (power transient and switch between electrolysis and fuel cell mode; ie. load cycling). Both cells and stacks have been optimized to increase performance (+20% compared to the one considered at the beginning of the project) in both SOFC and SOEC modes, To increase overall system efficiency by minimizing the energy losses in the Balance of Plant (BoP) components, specific reversible power electronics have been developed. With the support of modelling activities, the rSOC based system, made of 3 modules of 4 stacks each, coupled with a battery storage has been designed in a way to ensure achieving maximum efficiency, while still operating the system in safe and lifetime-optimal conditions. The site where the system will be installed has been prepared. Furthermore, techno-economic simulations have been carried out to evaluate market driven requirements on system sizing and costs to be competitive with other storage solutions.

Experimental

Cells developments

Different cell microstructures have been produced by Elcogen and compared with a reference structure. The reference structure was a fuel electrode supported cell, consisting of ~380 µm Ni/3YSZ support layer, a ~5 µm Ni fuel contact layer, a ~12 µm Ni/8YSZ fuel active layer, a ~7 µm thick LSC oxygen electrode, a ~2µm 8YSZ electrolyte and a CGO barrier layer of similar thickness. The modifications targeted the active fuel electrode, the fuel electrode support, the barrier layer and the oxygen electrode. Specifically, porosity, layer thicknesses and Ni/YSZ ratio were refined for the active and support fuel electrode structures whereas different thicknesses were produced for the oxygen electrode. The modification of the barrier layers aimed for the adjustment of the layer density, thickness and its manufacturing process. Additionally, a contact layer on oxygen electrode was added.

After testing those improvements at single cell level in terms of performance [START_REF] Ploner | [END_REF], two optimized cells were selected for further durability tests [START_REF] Hauch | 14th European SOFC&SOE Forum[END_REF] in rSOC conditions still at single cell (53x53 mm 2 ) level. Following the results obtained, and taking into account the ability of the cells to be produced in large amounts for the needs of the project and beyond, one of the two cells was considered as Generation 2 (G2) cell for the project.

Then full scale cells have been developed for integration into stacks for both lab tests at stack scale and for integration in the stacks to be embedded into the modules for in-field tests. Those cells are square, 120x120 mm, 100x100 mm (100 cm²) active area.

In parallel, to increase the power density further, enlarged cells have also been manufactured. The cell active area was nearly doubled passing from 100 cm² for reference cell to 196 cm², cells dimension being 160x160 mm 2 . For this work on enlarged cells; the reference microstructure was applied. Figure 1 illustrates the two sizes of cells considered for stacks in this project.

Stacks developments

Starting from the reference CEA 25-cell stack, comprising ~ 500 µm thick cells which was designed for a highly efficient SOEC operation [START_REF] Mougin | [END_REF], first an adaptation has been performed to integrate the reference ~ 400 µm thick aforementioned cells. Afterwards works have been performed to improve the fluidic distribution to the cells, especially on air side to cope with high flow rates as those needed for SOFC operation. Most of the work has been carried out to decrease the pressure drop on air side. Stack design modifications have been validated first at small scale, (5-cell short stack) and then at 25-cell scale considering reference cells [START_REF] Ploner | [END_REF].

A first 25-cell stack made of G2 cells was manufactured. In addition to the design optimizations described in [START_REF] Ploner | [END_REF] and to the integration of G2 cells, this stack integrates a manifold and an electrical insulation for an electrical connection in series of stacks into the modules for in-field test. After validation in lab of the improved behaviour of this G2-cell based stack, all stacks needed for the modules have been produced.

In parallel, 196 cm² active area cells have been integrated first in a small scale stack (5cell stack) and then a full-scale stack (25-cell stack) for performance validation.

For the purpose of the project, CEA produced five 5-cell stacks (four with 100 cm² cells and one with 196 cm² cells) and twenty-two 25-cell stacks (twenty-one with 100 cm² cells and one with 196 cm² cells). Out of those twenty-two stacks twelve stacks plus four spare stacks are dedicated to the integration into the modules for the in-field tests.

Cells and stacks testing

The overall protocols considered for cells and stacks testing at DTU and CEA are presented into REFLEX public deliverable D2.1 [START_REF]Deliverable D2.1, Cells and Stacks testing protocol[END_REF]. As an example, figure 2 presents the overall strategy for the performance and durability tests at stack level. For the rSOC cycling tests, two types of sequences were considered. Test #1 was conducted alternating from SOFC to SOEC mode by steps of ≈100 h, stack being fed with 50/50 H2O/H2 at total flow rate of 12 NmL/min/cm² on fuel side and air (clean and dry) on oxygen side. Tests were operated galvanostatic in SOFC and SOEC to reach respectively 0.8V (min) and 1.3 V (max) at test beginning, which corresponds to current density of 0.35 and -0.5 to -0.6 A/cm² at 700°C. Then test #2 was performed alternating from SOEC to SOFC by daily cycles of respectively 8h and 16h. Further details about the testing conditions can also be found in [START_REF] Ploner | [END_REF][START_REF] Hauch | 14th European SOFC&SOE Forum[END_REF][START_REF] Hauch | Load Cycling Tests of Reversible Solid Oxide Cells -Effects of Current Density, Steam Content and Utilization[END_REF][START_REF] Cubizolles | Stack optimization and testing for its integration in a rSOC-based renewable energy storage system[END_REF]. In addition, the stack performance has been tested in SOFC considering CH4 as a fuel (considered as representative of natural gas which will be considered for the real system in-field).

Modules design and modelling

The technology developed by Sylfen is called Smart Energy Hub (SEH). It is based on rSOC technology, that is to say able to operate either in SOEC (solid oxide electrolysis) mode or in SOFC (solid oxide fuel cell) mode in order to either store excess electricity to produce hydrogen, or when energy needs exceed local production, to produce electricity (and heat) again, from hydrogen or any other fuel locally available (natural gas, in particular). To be more effective, the rSOC system is completed with an electrochemical storage solution (Li-ion battery) allowing fast response to the electrical energy needs. The prototype designed in REFLEX intends to demonstrate the high power-to-power (P2P) round-trip efficiency of this technology (as compared to other H2 based solutions) and its flexibility and durability in dynamic operation (power transient and switch between electrolysis and fuel cell mode).

The system architecture (electric, thermal, electro-chemical, gas, …) has been defined by Sylfen. The product breakdown structure is given in Figure 3. -For the rSOC system: 3 modules of 4 stacks, with 3 operating modes: SOEC for H2 production, SOFC fed with H2, and SOFC fed with natural gas (NG) (table 1). For each mode 3 setpoints are defined, a minimum power Pmin, a nominal power Pmed and a maximum power Pmax. The maximum H2 production is 16 Nm 3 h -1 in SOEC mode. In SOFC mode, the maximum electrical power is up to 15 kW in SOFC, fed in H2 or NG. A standby mode and a safety mode are also defined. -For the batteries: a 50 kWh battery electrochemical storage system has been manufactured, including 7 Li-Al cells modules, with 4 cell packs per module and 144 unitary cells per pack. -Heat management is optimised in order to: 1) maximize system heat self-consumption, and 2) transfer a large amount of residual heat to the building hot water network.

Modelling tasks have been performed both to support the understanding and optimisation of the stack behaviour in terms of thermal management, to optimise the system design and to defined the best operation strategy. System modelling was performed using Simulink software utilizing the MATLAB-based ProSOFC component library developed by VTT. The system modeling activity was divided into 4 different phases, where each phase has objectives adapted to the current design level of the project: 1) concept level system design; 2) system modelling for design phase; 3) system modelling for operational phase; 4) interface to the end-user application. Results for the two first phases will be presented in next section.

The REFLEX Smart Energy Hub will be installed at ENVIPARK, in Torino, Italy. ENVIPARK is a Technology Park structured as a site of 10 building units that hosts offices and Laboratories operated by the various companies located in the Park. With real renewable energy sources upstream (a Hydroelectric power plant, 450 kWe power, 1700 MWh/y, and a Photovoltaic plant of 17 kW of peak power (13 MWh/y) integrated into the park's energy grid, and with real usages of electricity and heat downstream, the behaviour of the rSOC based Smart Energy Hub developed in REFLEX will be evaluated in real operating conditions.

Power electronics developments

Power electronics systems were manufactured by University of Seville and GPTech and developed according to the specification defined, taking into account two stacks in series in order to improve the DC/DC converter efficiency, and the battery modules, mainly the capacity and the peak charging and discharging powers.

The DC/DC converter system consists of six individual DC/DC modules of 20kW each, total system is 120kW power capable. Each one of the six individual DC/DC converters has the following structure (figure 4). The DC/DC converter is connected in the low voltage side to the rSOC system, that acts as a fuel cell (in SOFC mode) generating electrical power. Or in the other hand, it acts as an electrolyzer (in SOEC mode) generating hydrogen from electrical power. All modes are in the specifications table, corresponding to the setpoints targeted for the modules (Table 1). The DC/DC converters have been built and all the possible tests have been performed, confirming the correct operation in both SOFC and SOEC modes, according to the maximum power capabilities. Finally, AC/DC power converter and the associated control system have been tested with the DC/DC devices and the battery storage.

Techno-economic studies

Smart Energy Hubs integrate a rSOC technology, coupled with Li-ion batteries. Such a system has then to be seen in this context, as a power storage solution for buildings. It competes with two types of power storages: batteries and the combination of a low temperature water electrolyzer (alkaline or PEM) and a fuel cell.

The techno-economic study performed by ENGIE aims at evaluating the long-term market potential of this new solution developed by SYLFEN. A first phase was to study through simulations the best sizing of Smart Energy Hubs for several building markets in Europe. This meant here, defining the main parameters of the system (electrolyzer power input, fuel cell power output, hydrogen storage size, and battery capacity), assuming performance objectives for the technology (efficiencies). Based on those results, the requirements on CAPEX and OPEX (maintenance cost) have been defined for Smart Energy Hubs to be competitive with other storage technologies.

Four cases studies (building type) have been considered: offices, hotels, malls and aquatic centres. For each building type, one unique building has been considered, and studied over four different countries (France, Germany, Italy and United Kingdom). The maximum CAPEX and OPEX for Smart Energy Hubs have been calculated, based on simulations over one year. The requirements on CAPEX and OPEX have been defined for each reference solution that Smart Energy Hubs have to compete with that are: 1) Power and gas from the grid, and without PV nor power storage; 2) Power and gas from the grids, PV panels, without power storage, with or without export tariff; 3) Power and gas from the grid, PV panels, battery-based power storage; 4) Power and gas from the grids, PV panels, low-temperature hydrogen-based power storage.

In terms of OPEX, a five-years lifespan of rSOC stacks has been considered in this study, twice lower than those of low-temperature electrolyzer/fuel cell storage system and batteries. Stack replacement has been integrated here within annual OPEX.

Results

Cells and Stack Testing

Detailed results are presented in (5-6, 9-10). The purpose of the results presented here is to highlight the main achievements of the project as compared to the targets set. In terms of performance, REFLEX targets were to reach -1.2 A/cm² at 1.3V at 700°C in SOEC mode, and 0.6 A/cm² at 0.8V at 700°C in SOFC mode. Due to the developments performed by Elcogen for G2 cells, this goal has been achieved (figure 5), achieving high steam/fuel utilization (80% or even higher), highlighting a good electrode microstructure without or with minimum mass transport limitation. It corresponds roughly to 20% of improvement of current density at the thermoneutral voltage as compared to reference Elcogen cells tested in the same conditions. a) b) Figure 5. i-V curve of G2 cells recorded at 700°C; a) in SOEC, 90% H2O/10%H2; b) in SOFC in 100% H2 at fuel electrode, synthetic air at oxygen electrode; steam conversion in SOEC and fuel utilization in SOFC are reported on upper secondary x-axis

The benefit of G2 cells has also been confirmed at stack level (figure 6). Indeed, the maximum current density achievable without exceeding 1.4 V on the less performing cell is shifted from -1.0 A/cm² to -1.05 A/cm², and the scattering between the 25 cells of the stack is decreased, confirming the trend already reported at the scale of a 5-cell short stack [START_REF] Hauch | 14th European SOFC&SOE Forum[END_REF]. As compared to single cell the stack performance in similar conditions is a bit lower. It might be due to the stacking and the possible additional contact resistance which does not exist in single cell test setup, and/or to thermal effects, the exact temperature of the stack being more difficult to determine than for a single cell.

The operating window of the stack has been determined in SOEC and in SOFC fed with H2 and CH4 (representing in lab the NG case). Table 2 presents the values achieved, which will be the reference setpoints for the modules. They exceed the targets defined in the project, which were 50-100% in SOFC, and 70-100% in SOEC. Figure 6. Comparison of i-V curves recorded on 25-cell stacks with G2 cells and reference cells recorded at 700°C in SOEC, 90% H2O/10% H2 at fuel electrode, compressed air at oxygen electrode; steam conversion is reported on upper secondary x-axis 0.8 0.9

G2 cells

Ref cells

Nearly double active area cells have been successfully integrated into stacks, with the same level of performance achieved as with standard size cells [START_REF] Cubizolles | Stack optimization and testing for its integration in a rSOC-based renewable energy storage system[END_REF].

System modelling

Modelling activities have been performed to support the concept level system design. The system modelling has started on the basis of a rough description of the 3D structure of the planned system, considering all operating modes (SOEC, SOFC-H2 and SOFC-NG). Several options for the insulated hotbox have been tested in order to choose the best design to reach the overall performance. System modelling has shown that with a careful planning/construction it is possible to build an rSOC system capable of achieving the efficiency targets set. As an example, the figure 7 presents results for SOEC Pmax mode and shows that combinations of hotbox size and insulation thickness that meet SOEC efficiency requirements can be identified. This exercise has been repeated for all operating points and modes to select the best geometry possible for an operation of the SEH at the best performances on all operating windows. Subsequently, system modelling activities have supported the design phase to find a feasible system layout and dimensions including its BoP components to allow the system to meet its performance requirements. One key task was to determine the size of the operational windows for each mode. Figure 8 presents the operating window for the SOFC-H2 mode as an example, where heat management problems arise at high operating current densities. More detailed simulations related to the operational phase can be performed as the final 3D structure is designed and the final BoP components therein are experimentally characterized.

Figure 7. Combinations of hotbox size and insulation thickness that meet SOEC efficiency requirements (AC to HHV efficiency target of 80%), identified in green area. Figure 8. Operating windows for SOFC-H2 mode; black lines show width and location of the operational window (in current density) where the stack temperature can be controlled by the air flow rate for selected stack module heat loss estimate. The design value for heat loss from one four stack module is 600 W. The modelled design case is shown using red markers. The extrapolation trend is found by repeating the same analysis using different heat loss estimates for the same stack module Power electronics validation DC/DC converters were tested in both modes: SOFC and SOEC at their maximum power capabilities with the test bench shown in figure 9. First, tests were made to each DC/DC converter separately and the results obtained were compared. The results showed similar behavior for all converters. The best efficiency achieved in the laboratory for the DC-DC converter was 96% (figure 9) which is in line with the initial objective of the project (95%). Among the different types of buildings, the hotel case appears to be the most interesting for power storage, from a "technical" point of view. Indeed, the hotel case is the one which allows to consider a significant increase of building's autonomy thanks to power storage, as illustrated in Figure 10 for the Italian hotel case.

Figure 10. Origin of the power consumed by the hotel with PV panels covering 100% of roof surface, with and without storage, in Italy Another service that can be provided by the rSOC based power storage, is its ability to provide renewable heat to the building in addition to power. Indeed, only hydrogen-based storage is able to release heat that can be valorized within the building, through a water loop. For low-temperature electrolysis and fuel cell systems, the share of renewable heat remains however very limited (<1.5%). High-temperature rSOC system can contribute more significantly to the building's thermal needs, up to 7.5% for the considered case studies. However, a significant share of that heat provided by the storage comes from the operation of the system in natural gas CHP.

When considering reference solution 1), i.e. a building without storage and PV panels, the add-on installation of an rSOC based storage and PV panels allows to reach lower Total Cost of Ownership (TCO), as soon as rSOC module CAPEX is in the range of 2000 -5000 €/kW, depending on use case and country, except for France, where lower CAPEX are required due to lower power prices. The higher maximum CAPEX are achieved for the office, aquatic centre and shopping centre (figure 11).

Figure 11. Maximum CAPEX for the rSOC based storage solution to reach the same TCO of the same building without PV nor storage solution, depending on building's typology and country. These maximum CAPEX are expressed in €/kW of electrolysis, and refer only to the rSOC module, without H2 storage and batteries When considering reference solution 2), i.e. a building already equipped with a large PV surface, but without any power storage system, the add-on installation of a rSOC based storage system allows to reach lower TCO, as soon as rSOC module CAPEX is in the range of 1000 -1700 €/kW, depending on the use case, and mainly for countries as Germany, UK and Italy. While with reference 1, competitiveness is mainly achieved thanks to the installation of PV, with reference 2, which considers that PV panels are already installed, the addition of a power storage solution requires much quite low CAPEX, in order to achieve the same TCO.

Consideration of reference solution 3) is a way to compare directly the competitiveness of rSOC based solution with battery-based solution. For a building already equipped with PV panels, rSOC solution appears to be better from a TCO point of view than batteries, when rSOC CAPEX are lower than 2000 €/kW for the aquatic centre and the shopping centre, and in the range of 1500€/kW max for the hotel. The office use case shows more competitiveness for rSOC than batteries, with acceptable rSOC CAPEX in the range of 2000 to more than 6000 €/kW. Consideration of reference solution 4) is a way to compare directly the competitiveness of rSOC based solution with low temperature electrolysis and fuel cell-based storage solution. For a building already equipped with PV panels, rSOC solution appears to be better from a TCO point of view than low temperature electrolysis / fuel cell, when rSOC CAPEX are lower than 2000 €/kW for the aquatic centre, the shopping centre and the office, and in the range of 1000 -1500€/kW max for the hotel.

Sensitivity studies on the impact of power and gas prices showed that power price as a high and significant impact on rSOC system competitiveness whatever the reference solution is (figure 12). The higher the power price is, the more rSOC system is competitive with the reference solutions, allowing therefore higher CAPEX to equal or lower the TCO. The impact is the highest when compared to solutions without storage systems, as adding a storage when power price is high, allows to reduce the energy bill, thanks to a higher selfconsumption of PV. As rSOC system has also the ability to produce power from natural gas in CHP mode, its competitiveness is also increased with higher power prices versus battery and low-temperature electrolysis / fuel cell storage systems.

Figure 12. Sensitivity of maximum CAPEX for the rSOC solution with power price, depending on reference solution, for the shopping centre building in Germany. Natural gas price is set here at 0.038€/kWh Because of the CHP mode, higher gas prices, leading to reduced price difference with power prices, has a negative impact on rSOC competitiveness. Nevertheless, this impact remains quite limited.

Conclusion

In order to reach high power-to-power (P2P) round-trip efficiency, improved flexibility and durability in dynamic operation, the key components of the rSOC based Smart Energy Hub system have been optimised. Cells, stacks, and power electronics have been improved and targets set by the project achieved. The performance of cells have been increased by 20% as compared to reference cells, reaching -1.2 A/cm² at 1.3V in SOEC mode and 0.6 A/cm² at 0.8V in SOFC mode at 700°C at single cell level. The cells have been successfully integrated into stacks, which were tested in SOEC, SOFC in H2 and SOFC in CH4 to evaluate the operating window, which in terms of power covers the range 58-100, 23-100 and 13-100% respectively. Nearly double active area cells have been successfully integrated into stacks, with the same level of performance achieved as with standard size cells.

With the support of modelling activities, the rSOC based system, made of three modules of four stacks each, coupled with a battery storage has been designed. Several options for the insulated hotbox have been simulated in order to choose the best design to reach the overall performance. The optimum operating window has been defined. The site where the system will be installed has been prepared. Furthermore, techno-economic simulations have been carried out to evaluate market driven requirements on system sizing and costs to be competitive with other storage solutions.

Figure 1 .

 1 Figure 1. Picture of the 100 cm² and 196 cm² active area cells developed and delivered for integration into stacks in the project; left: 120x120 mm² cell; right: 160x160 mm² cell

  a): Stack testing flow chart; b): Mix gas composition set for i-V curves for 3 temperatures: 800, 750 and 700°C

Figure 3 .

 3 Figure 3. Product breakdown structure

Figure 4 .

 4 Figure 4. DC-DC converter electrical schematic

  9. a) Test bench used for lab validation of the DC/DC converters; b) Efficiency of the DC/DC converter recorded in SOEC mode as function of the current applied for three voltage values Techno-economic results

  (€/kW) to be competitive with same building without PV and storage FR IT UK DE

TABLE I .

 I DC-DC converter electrical operating specifications in agreement with the rSOC system setpoints.

	Operating mode		Current (A)	Voltage min (V)	Voltage max (V)	Power min (W)	Power max (W)
		Max	-130			7280	8840
	SOEC	Med	-80	56	76	4480	6080
		Min	-40			2240	3040
		Supermax	63			1512	3024
	SOFC-H2	Max Med	50 38	24	48	1200 912	2400 1824
		Min	24			576	1152
		Supermax	75			1800	3600
	SOFC-NG	Max Med	60 50	24	48	1440 1200	2880 2400
		Min	24			576	1152

TABLE II .

 II Operating window determined at the scale of a 25-cell stack made of G2 cells, given as a percentage of the max power

	Power (%)	SOEC	SOFC-H2	SOFC-CH4
	Pmin	58	23	13
	Pmed	80	66	75
	Pmax	100	100	100

Acknowledgments

This project has received funding from the Fuel Cells and Hydrogen 2 Joint Undertaking under grant agreement No 779577 (REFLEX project, http://www.reflexenergy.eu/). This Joint Undertaking receives support from the European Union's Horizon 2020 research and innovation programme and Hydrogen Europe and Hydrogen Europe Research.


Loïc Michel

Malek Ghanes

Franck Plestan

Yannick Aoustin

Jean-Pierre Barbot

Semi-Implicit Homogeneous Euler Differentiator for a Second-Order System: Validation on Real Data

Keywords: Discretization, Homogeneous differentiator, Attractor, Projectors, Experimental Data

In this paper, a semi-implicit Euler approximation scheme is proposed for a second-order homogeneous differentiator. Compared to an explicit Euler approximation, it is well-known that implicit Euler approximation scheme offers better performances like reducing high frequency oscillations. However, the implicit Euler approximation scheme works well only when dealing with classical sliding mode differentiator. In order to keep advantages of implicit Euler approximation, when this approximation is applied in case of homogeneous differentiators, a semi-implicit Euler approximation is proposed for a second-order system. Validation on real data is conducted to highlight the well-founded of the proposed differentiation strategy.

I. INTRODUCTION

This paper is dedicated to a real signal differentiation problem that is a crucial problem from a practical pointof-view; for example in robotics, a speed sensor is more expensive and fragile than an encoder sensor [START_REF] Spong | Robot dynamics and control[END_REF], [START_REF] Khalil | Modeling, identification and control of robots[END_REF]. Numerical differentiation is then preferred to design model based control. However, digital differentiation in the presence of Coulomb friction or noise can lead to instability or, at least, to high frequency phenomena of the closed loop system. This could be dramatic for example with a haptic control in a medical application. Improving the operation of differentiation is complex and therefore attracts the interest of a lot of scientists and a comparison of some of them has been performed in [START_REF] Mojallizadeh | Discrete-time differentiators in closedloop control systems: experiments on electropneumatic system and rotary inverted pendulum[END_REF]. The differentiation problem can be treated by two ways. The first way consists in using continuous-time differentiators [START_REF] Levant | Robust exact differentiation via sliding mode technique[END_REF] based on, for example, higher order sliding mode approaches and especially the well-known super-twisting algorithm [START_REF]Sliding order and sliding accuracy in sliding mode control[END_REF]. The second way is closer to real systems that are driven by discretized algorithms, and deals with differentiators based on numerical Euler approximation. In the literature, two approaches of Euler approximation can be used. The first one is a classical Euler explicit approximation approach, which has some drawbacks with respect to high frequency oscillations phenomena, see for example [START_REF] Mojallizadeh | Discrete-time differentiators in closedloop control systems: experiments on electropneumatic system and rotary inverted pendulum[END_REF]. Note that, in [START_REF] Livne | Proper discretization of homogeneous differentiators[END_REF], an explicit discrete differentiator including additional Taylor expansion corrective terms is proposed. {loic.michel, malek.ghanes,franck.plestan}@ec-nantes.fr Y. AOUSTIN is with Université de Nantes-LS2N, UMR 6004 CNRS, Nantes yannick.aoustin@univ-nantes.fr J.-P. BARBOT is also with QUARTZ EA 7393, ENSEA, Cergy-Pontoise, France barbot@ensea.fr Moreover, in [START_REF] Barbot | Discrete differentiators based on sliding modes[END_REF], the accuracy of such explicit discrete differentiator is studied. The second one, introduced in [START_REF] Acary | Implicit euler numerical scheme and chattering-free implementation of sliding mode systems[END_REF], is based on Euler implicit approximation approach and allows reducing significantly the oscillation phenomena whatever the actuator gain is [START_REF] Polyakov | Consistent Discretization of Finite-time and Fixed-time Stable Systems[END_REF]. In addition, a comparison of performances among 25 classical and advanced differentiation structures has been performed [START_REF] Mojallizadeh | Discrete-time differentiators in closedloop control systems: experiments on electropneumatic system and rotary inverted pendulum[END_REF], [START_REF] Mojallizadeh | Discrete-time differentiators: design and comparative analysis[END_REF] and highlights the fact that implicit based methods are more efficient and robust to the noise than the corresponding explicit strategies. Motivated by this ascertainment, this paper proposes a semiimplicit approximation, which is a combination between explicit and implicit techniques in order to be able to deal with homogeneous differentiators while keeping the benefits of Euler implicit approximation in the context of second order-system. The main contribution of the proposed differentiator is to combine explicit part with implicit one based on two projectors in order to reduce the effects of chattering phenomena as well as measurement noise and perturbations. To prove this claim, a theoretical study estimating the oscillation amplitude is firstly given in the context of explicit discretization in one dimensional case. To highlight the merit of the proposed strategy, experiments are conducted within an electrical RC circuit. This paper is outlined as follows. Section II describes the background of the homogeneous differentiation approach. Section III presents the main contribution and describes the semi-implicit homogeneous approach including a single projector, whereas the extension with two projectors is presented in Section IV. Section V describes a modification of the projectors to improve the robustness with respect to the noise. Both differentiation approaches are applied to experimental data in Section VI in which their benefits are discussed. Section VII concludes the paper and suggests some future works.

II. BACKGROUND ON HOMOGENEOUS DIFFERENTIATION

A. Basic analysis on the explicit Euler discretization

To introduce the homogeneous approach, consider the following continuous time system

ẋ = a (1)
where a is a constant perturbation. A continuous time estimator for (1) reads as (with b > 0)

ẋ = b x -x α (2)
where the notation

• α = | • | α sgn(•)
where the notation

• + = •(k + 1
) is adopted in the discretetime framework along the paper. For ease of reading, notation "e" is kept in the term of the right part of ( 4), but the term has also to be considered in the discrete-time framework; the equilibrium point of (4) being exactly the same than (3) i.e. e eq = a b 1 α . Remark that if α = 0, then e eq → 0 since b > a. In the sequel, two cases α ∈ (0, 1) and α > 1 are considered. The particular case α = 1 will be discussed at the end of this section.

CASE 1: 0 < α < 1.
It is well-known that e oscillates in the case of (4); considering the change of coordinate e = z + e eq yields to

z + = z + h(a -b z + e eq α ). (5)
A graph of z + versus z for ha = 10 -3 , hb = 1.7783 and α = 0.75 is displayed Fig. 1. From Fig. 1, quarters 1 and 3 appear to be unstable, whereas the quarters 2 and 4 are stable: as a consequence, contrarily to (3), the equilibrium point e eq is unstable for (4). Fig. 1 allows to define the domain Ω such that

Ω = {z, z + / z os1 ≤ z ≤ z os2 , z os1 ≤ z + ≤ z os2 }
where z os1 is the value of the abscissa of the point marked with a green diamond, and z os2 is the value of the abscissa of the point marked with a magenta diamond. Starting for example from the initial value denoted z 0 = 0.2 (marked with a black bullet), the system (5) converges to a limit cycle, whose behaviour oscillates between z os1 and z os2 such that

z os2 = -z os1 ≈ 0.6
as shown in Fig. 2. Now, with respect to Fig. 1, in order to evaluate the magnitude of the oscillation, it is necessary to evaluate a majorant of the non-zero coordinates of the intersection point, denoted z int , between z + -function (blue line -Fig. 1) and red lines of Fig. 1. The intersection points are defined such that

-z int = z int + h(a -b z int + e eq α )
which gives

2z int = hb z int + e eq α -ha. ( 6
)
Replacing e eq by a b 1 α , equation ( 6) becomes

2z int = hb z int + a b 1 α α -ha or 2z int = h( b 1 α z int + a 1 α α -a). (7)
Assuming that there exists w such that z int = w 1 α , equation ( 7) gives ( 8) and allows determining a majorant of z int depending on the following sign conditions:

2z int = h( bw 1 α + a 1 α α -a) (8)
1) if sgn( bw

1 α + a 1 α ) sgn(a) = 1, then two cases occur a) if | bw 1 α + a 1 α | α > |a| then, from
Jensen inequality dedicated to the convex functions [START_REF] Jensen | Sur les fonctions convexes et les inégalités entre les valeurs moyennes[END_REF], one has

| bw 1 α + a 1 α | α ≤ |bw| + |a|
and from ( 8), one gets

|2z int | ≤ h|bw| + |a| -|a| = h|bw| = hb|z int | α that gives the inequality |z int | = hb 2 1 1-α . b) if | (bw) 1 α + a 1 α | α ≤ |a|, it comes from (8) |2z int | ≤ |ha| then, for sgn( bw 1 α + a 1 α ) sgn(a) = 1, one has |z int | ≤ max hb 2 1 1-α , h|a| 2 . ( 9
)
2) if sgn( bw

1 α + a 1 α ) sgn(a) = -1, it implies |bw| 1 α > |a| 1 α . As a consequence, | bw 1 α + a 1 α | < |bw| 1 α
and |ha| < h|bw|. Then, equation [START_REF] Barbot | Discrete differentiators based on sliding modes[END_REF] gives

|2z int | < h|bw| + h|bw| ⇒ |z int | < hb|z int | α ⇒ |z int | < (hb) 1 1-α .
From ( 9) and the previous inequation, the worth case is

|z int | < (hb) 1 1-α . ( 10
)
CASE 2: α > 1. Figure 3 displays z + versus z for ha = 10 -3 , hb = 1.7783 and α = 2. It can be seen that the equilibrium point e eq is stable. Moreover, in Fig. 3, green and magenta diamonds highlight the points allowing to define the attractivity domain of e eq . Starting with an initial condition far from the equilibrium point, i.e. over the diamonds w.r.t. the origin of the (z, z + )plane, the explicit Euler discretization generates an unstable solution.

Remark 1: For the specific case α = 1, the curve is a line crossing the origin and either the line is included inside the quarters 1 or 3, and therefore (4) is unstable, or the line is included inside the quarters 2 and 4, and consequently (4) is globally stable. Moreover in ( 9) and [START_REF] Mojallizadeh | Discrete-time differentiators: design and comparative analysis[END_REF], it not possible to consider the case α = 1 because the limit goes to ∞.

B. Some recalls on two dimensional case Continuous-time system. Let p(t) be a bounded perturbation, which is considered unknown. The system under study consists of a double integrator defined as

           ẋ1 = x 2 ẋ2 = p(t) y = x 1 (11)
where x(t) ∈ R 2 is the state of the system, and y ∈ R is the output of the system.

Assumption 1: There exists p M > 0 such that |p(t)| < p M for all t > 0.

Homogeneous continuous-time differentiator. Given Assumption 1, from [START_REF] Perruquetti | Finite-time observers: application to secure communication[END_REF], [START_REF] Ghanes | A new varying gain exponent based differentiator/observer: an efficient balance between linear and sliding-mode algorithms[END_REF], a continuous-time homogeneous differentiator can be designed as follows

   ż1 = z 2 + λ 1 µ 1 α ż2 = λ 2 µ 2 1 2α-1 (12)
where 1 = x 1 -z 1 . The gains λ i > 0, i = {1, 2} and the parameter µ (sufficiently large to cancel the effect of p(t)) are tuned such that 1 and 2 = z 2 -x 2 converge towards zero. Then, the system (12) allows to get an estimation of ẋ1 thanks to z 2 .

Implicit Euler discrete-time system. Given Assumptions 1-2 on p(t) and h, the implicit Euler discretization of continuous-time system (11) reads as

   x + 1 = x 1 + h x + 2 = x 1 + h(x 2 + hp + ) x + 2 = x 2 + hp + . ( 13
)
Objective. The objective is now to propose an Euler discretization of the continuous-time homogeneous secondorder differentiator [START_REF] Perruquetti | Finite-time observers: application to secure communication[END_REF]. Two possibilities can be offered:

• Explicit homogeneous Euler differentiator This solution consists in proposing an explicit Euler approximation of the differentiator ( 12)

   x+ 1 = x1 + h (x 2 + λ 1 e 1 α ) x+ 2 = x2 + h λ 2 e 1 2 α-1 (14)
where e 1 = x 1 -x1 . However, this solution based on ( 14) is not attractive since it suffers of chattering phenomenon. • Implicit homogeneous Euler differentiator This solution reading as

   x+ 1 = x1 + h x+ 2 + λ 1 e + 1 α x+ 2 = x2 + h λ 2 e + 1 2 α-1 (15)
allows to propose an implicit Euler approximation [START_REF] Acary | Implicit euler numerical scheme and chattering-free implementation of sliding mode systems[END_REF] of [START_REF] Perruquetti | Finite-time observers: application to secure communication[END_REF] in the case of α = 0.5 (Levant differentiator).

It can be noticed that, from [START_REF] Porcelli | Bfo, a trainable derivative-free brute force optimizer for nonlinear bound-constrained optimization and equilibrium computations with continuous and discrete variables[END_REF], the observation error equation reads

e + 1 = e 1 + h e + 2 + λ 1 e + 1 α . ( 16
)
Remark that it is impossible to solve [START_REF] Qu | Fusion estimation for a class of multi-rate power systems with randomly occurring scada measurement delays[END_REF] under the constraint e + 1 = 0 for α = 0 due to the presence of |e + 1 | α associated to the correction term. Consequently, the implicit homogeneous Euler second-order differentiator [START_REF] Porcelli | Bfo, a trainable derivative-free brute force optimizer for nonlinear bound-constrained optimization and equilibrium computations with continuous and discrete variables[END_REF] does not work. That is the reason why, in the next section, a semi-implicit homogeneous Euler second-order differentiator is proposed, and its accuracy is evaluated considering experimental data.

III. SEMI-IMPLICIT HOMOGENEOUS EULER DIFFERENTIATOR WITH ONE PROJECTOR

A semi-implicit Euler discretization of the homogeneous differentiator [START_REF] Perruquetti | Finite-time observers: application to secure communication[END_REF] with a single projector (SIHD-1) is proposed and reads as

   x+ 1 = x1 + h x+ 2 + λ 1 µ|e 1 | α N 1 x+ 2 = x2 + E + 1 h λ 2 µ 2 |e 1 | 2α-1 N 1 (17)
where

N 1 :=      |e 1 | 1-α < λ 1 µh (e + 1 = 0) → N 1 = e1 1-α λ1hµ |e 1 | 1-α ≥ λ 1 µh (e + 1 = 0) → N 1 = sgn(e 1 ) (18)
is the projector that replaces the function sgn(e + 1 ) in [START_REF] Porcelli | Bfo, a trainable derivative-free brute force optimizer for nonlinear bound-constrained optimization and equilibrium computations with continuous and discrete variables[END_REF]; see [START_REF] Acary | Implicit euler numerical scheme and chattering-free implementation of sliding mode systems[END_REF] for detailed explanations.

Remark 2: As it is not possible to solve (16) because of the term |e + 1 | α sgn(e + 1 ), this latter is replaced in the proposed method by |e 1 | α N 1 , where N 1 formally replaces sgn(e + 1 ).

Moreover, E + 1 is defined as

E + 1 =    1 if e 1 ∈ SD 0 if e 1 / ∈ SD (19)
with SD defined as

SD = {e 1 / |e 1 | ≤ (λ 1 µh) 1 1-α }. ( 20
)
Remark 3: The terms |e 1 | α and |e 1 | 2α-1 are the explicit parts while the projector (the term N 1 ) refers to the implicit part. This projector N 1 cannot be evaluated for h = 0. Moreover, in (19), to take into account that, at the first time, e 1 ∈ SD, E 1 has no information with respect to e 2 (this will be explained later in (28)).

The differentiation error dynamics read as

   e + 1 = e 1 + h e + 2 -λ 1 µ|e 1 | α N 1 e + 2 = e 2 + h p + -E + 1 λ 2 µ 2 |e 1 | 2α-1 N 1 . ( 21
)
Before presenting the main results, the following assumption is needed.

Assumption 2:

1) there exists ẏM > 0 such that for all t > 0, | ẏ(t)| < ẏM ;

2) the perturbation p(t) is a constant parameter or slowly variable, that implies that for sufficient small h > 0, p + p. Remark 4: Assumption 2.1 is one of the conditions to stay in SD (see (26) and after); assumption 2.2 is an extra condition on the perturbation due to the Euler discretization and the fact that the variation of the perturbation must be negligible over a sampling period.

Theorem 1: Suppose that Assumptions 1 and 2 hold and he 2M ∈ SD. Then, for h > 0, there exist λ 1 > 0 and λ 2 > 0 such that the differentiation error dynamics (21) converge in finite time to

SD 1,2 = {e 1 , e 2 / e 1 ∈ SD 1 and e 2 ∈ SD 2 } with SD1 = {e 1 / |e 1 | ≤ max {he 2O1 , he 2O2 }} SD2 = {e 2 / |e 2 | ≤ max {e 2O1 , e 2O2 }} .
where

e 2O1 = p M λ 1 λ 2 µh α-1 - λ 2 µh α λ 1 (22)
and

e 2O2 = p M λ 1 λ 2 µh α-1 + λ 2 µh α 2λ 1 . ( 23
)
Proof: The proof is done in three steps: firstly the convergence in finite time of e 1 on SD is done using the stability analysis of (3). Afterwards, the equilibrium point of e 2 is calculated. Finally, the magnitude of the steady-state oscillations is deduced using the stability analysis of (3).

• Convergence of e 1 in (SD) . From ( 19)-( 20), e 1 ∈ SD if

|e 1 | 1-α < λ 1 µh ⇐⇒ |e 1 | < (λ 1 hµ) 1 1-α (24)
the conditions to reach SD are given in the sequel taking into account that e 1 = z + e eq . As it is expected at least that for z = z int , the solution is inside SD, and (24) reads

|z int | + |e eq | < (λ 1 hµ) 1 1-α (25)
1) if z has the same sign than e 2 , then from (5) including a = e 2 and b = λ 1 µ, considering (9) and the worth case |e 2 | = e 2M , the convergence in SD in finite time is given by the following inequality

e 2M λ 1 µ 1 α + hλ 1 µ 2 1 1-α < |λ 1 µh| 1 1-α (26)
where h > 0 and λ 1 µ > 0 are chosen to satisfy (26). Moreover, to stay in SD, e 2M used in (26) must be bounded. This condition is satisfied from assumption 2 and the fact that |x 2 | < ẏM then |e 2M | < 2 ẏM . 2) if z has not the same sign as e 2 , taking the same way than the previous case, considering [START_REF] Mojallizadeh | Discrete-time differentiators: design and comparative analysis[END_REF] and the worth case i.e |e 2 | = e 2M , it yields Remark that for e 2M = 0 the oscillation z is symmetrical around zero and the condition to reach SD is

|λ 1 µh| 1 1-α > max e 2M λ 1 µ 1 α -|hλ 1 µ| 1 1-α , e 2M λ 1 µ 1 α . ( 27
|λ 1 µh| 1 1-α > λ 1 µh 2 1 1-α
that is always verified for h, λ 1 and µ strictly positive. When e 1 ∈ SD, its implies that

E + 1 = 1 and N 1 = e 1 1-α λ 1 hµ ,
then e + 1 verifies the following equation

e + 1 = he + 2 (28)
and e 1 stays on SD if |e + 2 | is smaller than e 2M , which will be proved in the sequel. Then, from (28), on SD, the second row of ( 21) becomes

e + 2 = e 2 + h p + - λ 2 µ λ 1 h he 2 α . ( 29
)
Remark 5: Equation ( 29) is the explicit Euler discretization of the following continuous dynamics

ė2 = p + - λ 2 µh α-1 λ 1 e 2 α ( 30
)
and it is well known that explicit Euler discretization of (30) oscillates.

• Equilibrium Point of (29). The equilibrium point of ( 29) is the solution of

p + - λ 2 µ λ 1 h he + 2 α = 0. ( 31
)
This gives the following equilibrium point e 2,eq for e 2 e 2,eq = λ 1 λ 2 µ p

+ 1 α h 1-α α . (32)
• Steady-state oscillations of e 2 and e 1 . Remark that (29) is exactly ( 4) with b = λ2µh α-1 λ1 and a = p + . Then setting e 2 = z + e 2eq , equation ( 5) is obtained with e eq = e 2eq . From the previous analysis, e 2 converges in the domain

SD2 = {e 2 / |e 2 | ≤ max {e 2O1 , e 2O2 , }} , (33)
where e 2O1 and e 2O1 are respectively defined in ( 22)-(23).

Remark that the absolute values are removed because all the terms are positive, and e 2M is chosen larger than e 2O1 and e 2O2 . Finally, from (28), the oscillation within e 1 is confined inside the SD1 domain

SD1 = {e 1 / |e 1 | ≤ max {he 2O1 , he 2O2 }} . (34)

IV. SEMI-IMPLICIT HOMOGENEOUS EULER DIFFERENTIATOR WITH TWO PROJECTORS

In order to remove the oscillation on e 2 and consequently on e 1 , the following differentiator with two projectors (SIHD-2) is proposed

     x+ 1 = x1 + h x+ 2 + λ 1 µ|e 1 | α N 1 x+ 2 = x2 + E + 1 h λ 2 µ 2 |e 1 | 2 α-1 N 2 (35)
where x+ 2 is computed in the first and second rows. The associated projectors are defined by

N 1 :=      |e 1 | 1-α < λ 1 µh (e + 1 = 0) → N 1 = e1 1-α λ1hµ |e 1 | 1-α ≥ λ 1 µh (e + 1 = 0) → N 1 = sgn(e 1 ) (36)
and, given that, on SD, e 1 = he 2

N 2 :=      |e 1 | 2-2α < λ 2 µ 2 h 2 → N 2 = e1 2-2α λ2h 2 µ 2 |e 1 | 2-2α ≥ λ 2 µ 2 h 2 → N 2 = sgn(e 1 ) (37)
with E 1 defined in (19). The estimation error dynamics read as

     e + 1 = e 1 + h e + 2 + λ 1 µ|e 1 | α N 1 e + 2 = e 2 + h p + -E + 1 λ 2 µ 2 |e 1 | 2 α-1 N 2 . (38)
Theorem 2: Suppose that assumptions 1-2 hold. Then for h > 0, there exist λ 1 > 0 and λ 2 > 0 such that the differentiation error dynamics (38) converge in finite time to

SD 1,2 = {e 1 , e 2 / e 1 ∈ SD 1 and e 2 ∈ SD 2 } with SD 1 = e 1 / |e 1 | ≤ h 2 p M , SD 2 = {e 2 / |e 2 | ≤ hp M } .
Proof: The beginning of the proof is similar to the proof of Theorem 1 concerning the convergence on SD. After this convergence, E 1 = 1, and since e 1 = he 2 , the dynamic of e + 2 becomes

e + 2 = e 2 + h p + -λ 2 µ 2 h 2α-1 |e 2 | 2 α-1 N 2 .
Choosing λ 2 sufficiently large after some iterations, it reads |e 2 | 2-2α < λ 2 µ 2 h 2α . Then, substituting again e 1 = he 2 , the projector N 2 can be read as V. θ PROJECTORS To improve the robustness with respect to the noise, a modification could be made for the projectors N 1 , and N 2 as follows

N 2 = e2 2-2α λ2h
N 1 :=      (1 -θ)|e 1 | 1-α < λ 1 µh → N 1 = (1-θ) e1 1-α λ1hµ (1 -θ)|e 1 | 1-α ≥ λ 1 µh → N 1 = sgn(e 1 )
and

N 2 :=      |(1 -θ)e 1 | 2-2α < λ 2 µ 2 h 2 → N 2 = (1-θ) e1 2-2α λ2h 2 µ 2 |(1 -θ)e 1 | 2-2α ≥ λ 2 µ 2 h 2 → N 2 = sgn(e 1 )
with θ ∈ [0, 1). Indeed, to damp the influence of the noise on the sliding surface, instead of requesting e + 1 = 0, it is only requested e + 1 = (1 -θ)e 1 . From the stability perspective, considering θ equal to zero, finite-time convergence is ensured and when θ is different from zero, the differentiator is less sensitive to noise, providing only an asymptotic convergence, since the influence of noise on the sliding surface is multiplied by (1 -θ).

VI. EXPERIMENTAL APPLICATION

In this section, a comparative study illustrates the properties of the proposed versions of the semi-implicit homogeneous differentiation technique. Similarly to the previous work of the authors [START_REF] Michel | An experimental investigation of discretized homogeneous differentiators: pneumatic actuator case[END_REF] where some homogeneous differentiation algorithms have been compared in the framework of velocity/acceleration estimation regarding the position of a pneumatic actuator, this study deals with the estimation of the current in the case of an electronics RC low pass circuit. It aims to compare the performances of the proposed differentiation solutions with a standard first order Butterworth filter 1 . This comparison is mainly based on the analysis of the Sum of Square Error (SSE)2 applied to the error between each differentiation version and the measured current.

The objective is to differentiate the output of a RC series circuit (R = 100 Ω and C = 100 µF), this signal being a sine function of frequency ω = 64.71 rad/s (f = 10.3 Hz). The measured signal to differentiate -see Fig. Two sampling periods h 1 = 2 10 -4 s and h 2 = 2 10 -3 s are considered in the sequel to highlight the properties of the SIHD method at lower sampling.

The Butterworth filter is set at the cutoff frequency of 100 Hz; the parameters of the SIHD methods λ 1 = λ 2 = 10 6 are chosen to have the linear part stable, and the θ 1 , θ 2 , α and µ parameters have been obtained thanks to an optimization procedure in order to improve the SSE performance index [START_REF] Porcelli | Bfo, a trainable derivative-free brute force optimizer for nonlinear bound-constrained optimization and equilibrium computations with continuous and discrete variables[END_REF]. These choices are summarized in the following table. 39)) and the outputs of the three differentiators (Butterworth, and the both proposed ones) for the two sampling periods h 1 and h 2 . It appears that, globally, these three solutions allow to get an image of the differentiation, but with a huge difference of performances. Furthermore, the value of the sampling period has an impact on the quality of the differentiation. SIHD-1 and SIHD-2 according to the Butterworth filter. For both samplings h 1 and h 2 , the semi-implicit homogeneous differentiators show better accuracy regarding the SSE index than the classical Butterworth filter. Since the two-projectorbased version (SIHD-2) includes a dedicated projector to drive the tracking error on e 2 and suppress the corresponding oscillations, it provides slightly better performances than the corresponding SIHD-1 method. Moreover, the SIHD-2 method is easier to extend to higher order. Nevertheless in applications, performances and the feasibility of both methods are quite similar (see [START_REF] Michel | An experimental investigation of discretized homogeneous differentiators: pneumatic actuator case[END_REF]).

h [s] method θ 1 θ 2 α µ SSE 2.10 -4 Butterw. - - - - 2

VII. CONCLUSIONS AND FUTURE WORKS

This paper proposes a semi-implicit Euler approximation of an homogeneous differentiator for a second-order system. The proposed scheme allows to apply an implicit Euler approximation (combined with explicit one) when homogeneous differentiators are considered instead of classical sliding mode differentiators. Two versions of differentiators based on one or two projectors are proposed to reduce the effect of chattering phenomena, measurement noise and perturbation. These two differentiators and a discrete Butterworth filter (used for comparison) have been implemented and tested on an experimental electronics RC filter. Better performances have been obtained within the proposed two differentiators compared to the discrete Butterworth filter in term of tracking and noise rejection. Future investigations include the generalization of the proposed approaches to higher order differentiation as well as extension to variable exponent differentiators to deal with adaptive noise rejection [START_REF] Ghanes | A new varying gain exponent based differentiator/observer: an efficient balance between linear and sliding-mode algorithms[END_REF]. In addition, more comparisons and investigations will be conducted with respect to recent implicit differentiation techniques. Moreover, it will be very interesting to adapt the Taylor expansion corrective term proposed in [START_REF] Livne | Proper discretization of homogeneous differentiators[END_REF] to the proposed implicit methods and also to consider not only a regular sampling but also the self and event triggering case [START_REF] Barbot | Discrete differentiators based on sliding modes[END_REF], [START_REF] Qu | Fusion estimation for a class of multi-rate power systems with randomly occurring scada measurement delays[END_REF].

1 α

 1 is adopted all along the paper. So, the dynamics of the estimation error e = x-x read as ė = a -b e α (3) Setting α ∈ (0, 1), and supposing that b is sufficiently large versus a, i.e. b |a|, e converges to e eq = a b in a finite time t f > 0. The explicit Euler discretization of (3) gives e + = e + h(a -b e α ),

Fig. 1 .

 1 Fig. 1. z + versus z (blue line) for α = 0.75.

Fig. 2 .

 2 Fig. 2. z versus t (sec) for z(0) = 1.

Fig. 3 .

 3 Fig. 3. z + versus z (blue line) for α = 2.

  ) It is always possible to find λ 1 µ such that (27) is verified.

2α µ 2

 2 and consequently e + 2 = hp + and finally |e 2 | ≤ hp M and |e 1 | ≤ h 2 p M .

Remark 6 :

 6 4 has been obtained by scoping the voltage v c across the capacitor C.Moreover, the voltage v r = R i c is accordingly scoped across the terminals of the resistor R, and constitutes the reference measured differentiated signal since By a practical point-of-view, the interest of the differentiator is to avoid the use a current sensor given that i c = C vc .

Fig. 4 .

 4 Fig. 4. Measured vc versus time (s) for h 2 = 2 10 -4 s.

Figures 5 -

 5 Figures 5-6 display the real differentiated signal (obtained from measurement -see (39)) and the outputs of the three differentiators (Butterworth, and the both proposed ones) for the two sampling periods h 1 and h 2 . It appears that, globally, these three solutions allow to get an image of the differentiation, but with a huge difference of performances. Furthermore, the value of the sampling period has an impact on the quality of the differentiation.

Fig. 5 .

 5 Fig. 5. Differentiated signals (measured, Butterworth, proposed methods) versus time (s) for h 1 = 2 10 -4 s.

Figures 7

 7 Figures 7 and 8 compare the SSE performances 3 of the

Fig. 6 .Fig. 7 .

 67 Fig. 6. Differentiated signals (measured, Butterworth, proposed methods) versus time (s) for h 2 = 2 10 -3 .

Fig. 8 .

 8 Fig. 8. Evaluation of the normalized index SSE h 2 = 2 10 -3 s.

The butter function of Matlab R is used to perform a discrete first order high pass filtering.

The SSE is given by SSE(•) = 1 n n l=1 • 2l with n the data size; the second half of the signal (in the steady state) has been considered for the SSE evaluation.

The normalization of SSE in Figures7 and 8is made with respect to results obtained with the Butterworth filter (always equal to one).

VIII. ACKNOWLEDGMENTS

This work is supported by Agence Nationale de la Recherche (ANR) with project DIGITSLID ANR-18-CE40-0008-01.


Pulse-level Noise Mitigation on Quantum Applications Siyuan Niu * and Aida Todri-Sanial † LIRMM, University of Montpellier, France (Dated: April 4, 2022)

Currently available quantum computers are prone to errors. Circuit optimization and error mitigation methods are needed to design quantum circuits to achieve better fidelity when executed on NISQ hardware. Dynamical decoupling (DD) is generally used to suppress the decoherence error and different DD strategies have been proposed. Moreover, the circuit fidelity can be improved by pulselevel optimization, such as creating hardware-native pulse-efficient gates. This paper implements all the popular DD sequences and evaluates their performances on IBM quantum chips with different characteristics for various well-known quantum applications. Also, we investigate combining DD with pulse-level optimization method and apply them to QAOA to solve Max-Cut problem. Based on the experimental results, we found that DD can be a benefit for only certain types of quantum algorithms, while the combination of DD and pulse-level optimization methods always has a positive impact. Finally, we provide several guidelines for users to learn how to use these pulse-level noise mitigation methods to build circuits for quantum applications with high fidelity.

I. INTRODUCTION

Quantum computing is rapidly growing in recent years and various technologies have been developed for different quantum platforms, the leading candidates being superconducting and trapped-ion devices. Several companies such as IBM, Rigetti, and IonQ provide publicly available cloud-based services which allow users to access their platforms remotely. But, today's quantum computers are still prone to errors caused by either unavoidable interactions with the environment or imperfect quantum controls. They are qualified as Noisy Intermediate-Scale Quantum (NISQ) computers [START_REF] Preskill | Quantum computing in the nisq era and beyond[END_REF]. The largest quantum chip to date has 127 qubits released by IBM [START_REF]IBM quantum unveiled eagle, a 127-qubit quantum processor[END_REF]. Quantum error correction (QEC) has been proposed to eliminate the noise impact and help achieve a fault-tolerant quantum device [START_REF] Calderbank | Quantum error correction via codes over gf (4)[END_REF][START_REF]Quantum error correction for quantum memories[END_REF]. However, the implementation of quantum error correction codes requires a large number of ancilla qubits, which is not feasible on current hardware. Therefore, alternative approaches are needed to tackle the noise issue.

In order to reduce the noise impact and improve the quantum circuit fidelity, several quantum software were designed and they have made a great contribution to different circuit design processes, such as circuit synthesis [START_REF] Matteo | Parallelizing quantum circuit synthesis[END_REF][START_REF] Zhang | Optimal quantum circuit synthesis from controlled-unitary gates[END_REF], or qubit mapping [START_REF] Niu | A hardware-aware heuristic for the qubit mapping problem in the nisq era[END_REF][START_REF] Li | Tackling the qubit mapping problem for nisq-era quantum devices[END_REF]. In addition, quantum error mitigation (QEM) was introduced for error suppression on NISQ devices [START_REF] Endo | Practical quantum error mitigation for near-future applications[END_REF]. There are different QEM techniques, such as readout error mitigation [START_REF] Bravyi | Mitigating measurement errors in multiqubit experiments[END_REF], dynamical decoupling (DD) [START_REF] Souza | Robust dynamical decoupling[END_REF], crosstalk mitigation [START_REF] Murali | Software mitigation of crosstalk on noisy intermediate-scale quantum computers[END_REF][START_REF] Niu | Analyzing crosstalk error in the nisq era[END_REF], zeronoise extrapolation (ZNE) [START_REF] Li | Efficient variational quantum simulator incorporating active error minimization[END_REF], etc. Most methods require supplementary circuit executions to build the error map for mitigation. Whereas DD is one of the simplest strategies, which aims at mitigating the decoherence error (also called idle error) without any circuit overhead, and is the focus of this paper. The thrust of DD is to insert periodically a series of pulses to the idle qubits and return the qubits to their original states. There is a high probability to have idle qubits during execution due to the variation of gate latencies and limited parallelism caused by the anti-commutative gates. It has been shown in [START_REF] Das | Adapt: Mitigating idling errors in qubits via adaptive dynamical decoupling[END_REF] that the idle qubits can be almost 10 times more subject to errors if adjacent two-qubit operations are executing at the same time on IBM superconducting device. DD plays an important role in reducing the idle error and has been used in quantum volume experiments [START_REF] Jurcevic | Demonstration of quantum volume 64 on a superconducting quantum computing system[END_REF], noise spectrum characterization [START_REF] Sung | Non-gaussian noise spectroscopy with a superconducting qubit sensor[END_REF], and decoherenceprotected quantum gate implementation [START_REF] Van Der Sar | Decoherence-protected quantum gates for a hybrid solid-state spin register[END_REF], etc. Moreover, there are different DD strategies, such as Hahn echo [START_REF] Solomon | Multiple echoes in solids[END_REF], CPMG [START_REF] Meiboom | Modified spin-echo method for measuring nuclear relaxation times[END_REF], XY4 [START_REF] Viola | Dynamical decoupling of open quantum systems[END_REF], robust KDD [START_REF] Souza | Robust dynamical decoupling for quantum computing and quantum memory[END_REF], etc. They have diverse impacts on decoherence error suppression for different quantum devices. It has been shown in [START_REF] Pokharel | Demonstration of fidelity improvement using dynamical decoupling with superconducting qubits[END_REF] that DD is able to extend the lifetime of one-qubit states as well as entangled two-qubit states for IBM and Rigetti devices using XY4 sequence. But this sequence is shown to be vulnerable to experimental imperfections, while the use of robust DD was demonstrated to be capable of correcting the pulse errors on Rigetti devices [START_REF] Souza | Process tomography of robust dynamical decoupling with superconducting qubits[END_REF].

Recently, IBM released Qiskit Pulse [START_REF] Mckay | Qiskit backend specifications for openqasm and openpulse experiments[END_REF][START_REF] Alexander | Qiskit pulse: programming quantum computers through the cloud with pulses[END_REF] allowing users to design and customize the gate pulse implementations. Some works have already attempted to optimize the pulse controls and reduce the pulse durations with the help of Qiskit Pulse [START_REF] Shi | Optimized compilation of aggregated instructions for realistic quantum computers[END_REF][START_REF] Khaneja | Optimal control of coupled spin dynamics: design of nmr pulse sequences by gradient ascent algorithms[END_REF][START_REF] Egger | Adaptive hybrid optimal quantum control for imprecisely characterized systems[END_REF]. However, such method often needs additional calibrations which is time-consuming and requires a deep familiarity with quantum control. In [START_REF] Stenger | Simulating the dynamics of braiding of majorana zero modes using an ibm quantum computer[END_REF][START_REF] Earnest | Pulseefficient circuit transpilation for quantum applications on cross-resonance-based hardware[END_REF], the authors proposed a new technique to create more hardware-native pulse-efficient gates, which improves the gate fidelity without the overhead of extra calibrations.

So far, the aforementioned techniques were only tested separately for limited benchmarks and quantum hardware. Several questions still remain unclear: [START_REF] Preskill | Quantum computing in the nisq era and beyond[END_REF] What are the impacts of different DD sequences on specific quantum algorithms? (2) For a certain benchmark, does the impact of different DD sequences vary across differ-ent quantum chips? (3) Will the combination of DD and pulse-level optimization methods further improve the circuit fidelity? In our work, we address these questions and our main contributions can be listed as follows:

• We explore the performance of different DD sequences on various quantum applications and evaluate the experiments on several IBM devices with different qubit numbers and quantum volumes. To the best of our knowledge, this is the first attempt to illustrate the behavior of applying robust KDD sequence to IBM quantum devices.

• We combine the DD technique with pulse-efficient gate optimization method to demonstrate their benefits on Quantum Approximate Optimization Algorithm (QAOA) to solve Max-Cut problem.

• Based on the experimental results, we provide guidelines and insights for users to realize application-oriented pulse-level error mitigation.

II. BACKGROUND A. Dynamical decoupling

Dynamical decoupling is wildly used in suppressing the decoherence error by reducing the interaction between the system and the environment. Considering a systembath Hamiltonian H shown in (1), H s and H B being the Hamiltonian of the system and the bath respectively, and the interacting term being H SB (see (2)), where σ α i is the Pauli matrix acting on qubit i, B α i is the operator of the environment, and α ∈ {x, y, z}. DD aims at reducing the impact of system-environment interaction and various DD protocols have been developed to improve the performance of quantum computers. Here, we review the main strategies of DD implementations and a summary of these DD sequences is shown in Table I.

H = H s + H B + H SB (1)
H SB = N i=1 α∈{x,y,z} σ α i ⊗ B α i (2)
• Hahn echo. The spin-echo sequence [START_REF] Solomon | Multiple echoes in solids[END_REF] is used to reduce the inhomogeneous effects from the environmental magnetic field. It applies a π pulse to the spin system to inverse the spins after a period of time t and let the system be refocused during the same duration. At that time, T * 2 effects are removed and we can obtain a T 2 echo.

(π) π 6 +φ -(π) φ -(π) π 2 +φ -(π) φ -(π) π 6 +φ
where (θ) φ = e • CP and CPMG. Carr and Purcell proposed a series of π pulses separated by a constant interval, known as the CP sequence, to further reduce the effect of self-diffusion in the inhomogeneous field [START_REF] Carr | Effects of diffusion on free precession in nuclear magnetic resonance experiments[END_REF]. But, extra pulses can introduce more errors and destroy the state of the system. Therefore, Meiboom and Gill improved the CP pulse by developing the CPMG sequence [START_REF] Meiboom | Modified spin-echo method for measuring nuclear relaxation times[END_REF], retaining the CP pulse but introducing an additional phase shift in the first pulse to compensate the pulse errors.

• UDD. Hahn echo and CP/CPMG are equidistant pulses. While UDD aims at optimizing the π sequence based on CP/CPMG by varying the intervals between each pulse [START_REF] Uhrig | Keeping a quantum bit alive by optimized π-pulse sequences[END_REF]. It is proven to be the optimal pulse to suppress low frequency noise and is insensitive to thermal fluctuations. It outperforms the equidistant DD sequences especially for systems whose spectral densities have high frequencies with sharp cutoff. For other more general cases, equidistant DD sequences were shown to perform better [START_REF] Ajoy | Optimal pulse spacing for dynamical decoupling in the presence of a purely dephasing spin bath[END_REF].

• XY4. All the DD sequences above only rotate around one single axis. They are exclusively useful when the system-environment interaction is orthogonal to the rotation axis. XY4 is the simplest DD sequence to generally suppress system-environment interaction along three directions [START_REF] Viola | Dynamical decoupling of open quantum systems[END_REF]. It inserts alternatively π rotations around X and Y axes independently of the initial state.

• XY8 and XY16. DD can sometimes detriment the fidelity due to the accumulation of errors caused by pulse imperfections. One approach to reducing the errors is to combine the basic sequence with its inverse for self-correction so that the unwanted terms can be canceled [START_REF] Gullion | New, compensated carr-purcell sequences[END_REF]. XY4 is usually chosen as the basic cycle. XY8 is composed of XY4 se-quence and its inverse while XY16 contains XY8 and its inverse.

• KDD. The other approach to avoid pulse imperfections is to replace each pulse of a DD sequence by a robust composite pulse, which is designed to generate ideal rotations even if there exist pulse imperfections [START_REF] Grant | Encyclopedia of nuclear magnetic resonance[END_REF]. KDD constructs a DD block using a 5-pulse composite π pulse, and combine two of such DD block with and without a phase shift π/2 as (KDD φ , KDD φ+π/2 ) 2 [START_REF] Souza | Robust dynamical decoupling for quantum computing and quantum memory[END_REF]. It is composed of 20 pulses in total.

DD has been demonstrated to have the capability to mitigate the decoherence errors on IBM and Rigetti platforms using XY4 sequences [START_REF] Pokharel | Demonstration of fidelity improvement using dynamical decoupling with superconducting qubits[END_REF]. The experiment is performed as follow: it prepares different initial states by varying the angles of the rotation gates, inserts XY4 sequences, and compares the difference between input and output states. The impacts of different DD sequences on Rigetti device are reported in [START_REF] Souza | Process tomography of robust dynamical decoupling with superconducting qubits[END_REF]. A large number of DD sequences are evaluated, including XY4, XY8, KDD, etc., and quantum process tomography is used to characterize the evolution. KDD is shown to be the most robust pulse sequence against pulse imperfections. Both papers evaluate DD sequences on simple benchmarks instead of real quantum applications, and there is no experiment to illustrate applying KDD to IBM quantum devices. Moreover, the experiments in [START_REF] Das | Adapt: Mitigating idling errors in qubits via adaptive dynamical decoupling[END_REF] show that the naive implementation of DD (inserting DD sequences to all the idle qubits when it is possible) can not always improve the circuit fidelity. Therefore, Das et al. proposed an Adaptive Dynamical Decoupling framework to estimate the DD impact for each circuit and adjust DD sequence to ensure it improves the circuit fidelity [START_REF] Das | Adapt: Mitigating idling errors in qubits via adaptive dynamical decoupling[END_REF]. This method achieves fidelity improvement but introduces a large overhead of DD impact characterization for a given application. Also, Ravi et al. proposed VAQEM [START_REF] Ravi | Vaqem: A variational approach to quantum error mitigation[END_REF], an approach that dynamically inserts DD sequences for variational algorithm with the overhead of tuning DD sequence. In our paper, instead of carefully adjusting DD sequences for better error mitigation with extra circuit overhead, we exhaustively evaluate DD sequences on extensively used quantum applications, to give a highlevel idea about if general DD techniques can really be beneficial for these applications. It is also the first attempt to implement KDD on IBM devices to check its performance. The experiments are performed on various IBM devices with different qubit numbers and quantum volumes to provide general insights about applicationoriented DD noise mitigation on IBM quantum devices.

B. Pulse-efficient technique

CN OT gate is the only two-qubit operation included in the basis gates for IBM quantum devices and its calibrations are provided through the IBM quantum dashboard. It is implemented by a R ZX ( π 2 ) gate along with some single-qubit gates [START_REF] Chow | Implementing a strand of a scalable fault-tolerant quantum computing fabric[END_REF] on the device. R ZX gate is realized by the echoed cross-resonance gate [START_REF] Sheldon | Procedure for systematically tuning up crosstalk in the cross-resonance gate[END_REF], which is specific for the IBM fixed-frequency superconducting transmon qubit device due to its low overhead and high fidelity. When executing a circuit on IBM quantum hardware, every multi-qubit operation needs to be transpiled to the CN OT basis which is not flexible and less efficient. Therefore, [START_REF] Stenger | Simulating the dynamics of braiding of majorana zero modes using an ibm quantum computer[END_REF][START_REF] Earnest | Pulseefficient circuit transpilation for quantum applications on cross-resonance-based hardware[END_REF] proposed a pulse-efficient circuit transpilation framework and has shown to achieve higher fidelity than CN OT basis transpilation for certain benchmarks. Instead of only using a R ZX (θ) with a fixed angle θ = π 2 for CN OT gate, a flexible echoed R ZX gate is implemented with arbitrary angle and R ZX (θ) = XR ZX (-θ

2 )XR ZX ( θ 2 ) to enable transpilations to R ZX basis. This method does not require any additional pulse calibration since the calibration of R ZX gate can be easily calculated from the CN OT calibrations. An example of transpiling R ZZ gate to two different basis gate sets is shown in Fig. 1. We first transpile the R ZZ gate to the IBM CN OT -basis gate sets (Fig. 1(b)) and then transpile again the original circuit to R ZX -basis gate sets (Fig. 1(c)). The duration of the second circuit is reduced by 51.1% compared to the first gate sets. This technique exposes the echo of the cross-resonance gates which enables at most one single-qubit gate between each non-echoed R ZX gate so that it can shorten the total circuit duration.

Quantum algorithms requiring a lot of R ZZ gate implementations can be particular benefited from this pulseefficient transpilation framework, such as QAOA [START_REF] Farhi | A quantum approximate optimization algorithm[END_REF], quantum hamiltonian simulation [START_REF] Papageorgiou | On the efficiency of quantum algorithms for hamiltonian simulation[END_REF], etc.

III. METHODS

We construct DD sequences according to Table I, and n is the number of repetition time for the basic DD pulse. For Hahn echo, we insert one X or Y gate during the idle time. We set n = 2 for CP/CPMG, CP (resp.CPMG) be-

q 0 : H • H q 1 : H • H q 2 : H • H q 3 : X H (a) BV algorithm q 0 : H X • X H • H q 1 : H X • X H • H q 2 : H X • X H • H q 3 : H • H • H (b) HS algorithm q 0 : H • q 1 : H • • q 2 : H • • q 3 : H • • q 4 : H • (c) Graph State q 0 : • • • P ( π 2 ) H × q 1 : • • P ( π 2 ) H P ( π 4 ) • × q 2 : • P ( π 2 ) H P ( π 4 ) • P ( π 8 ) • × q 3 : H • • • × (d) QFT
FIG. 2: Examples of different quantum circuits. Note that, the Graph State circuit is constructed linearly based on the target hardware topology.

ing implemented as the sequence X -X (resp. Y -Y ). For XY4, XY8, and XY16, n is set to 1 so that the sequence of pulses is inserted once. All the aforementioned DD sequences are equidistant and symmetrical: there is a delay of τ /2 at the beginning and the end of the idle time, and a delay of τ between each DD pulse. Whereas UDD is composed of non-equidistant pulses, it acts differently from CPMG for all n > 2. As there is no clear conclusion about which repetition number n we should pick to have the best performance of UDD on IBM quantum devices, we set n = 8 as demonstrated in Qiskit tutorial. The repetition gate is set to X and Y , and is marked as UDD X and UDD Y in the results. For KDD, (θ) φ represents a rotation of θ around the axis defined by φ.

Since θ φ = R Z (φ)R X (-θ)R Z (-φ), we insert 60 gates in total, where delay is only introduced between each θ φ , which means there is no space between the three gates constructing θ φ . Note that, since R Z gate is free on IBM quantum devices, KDD is translated to "only" 20 pulses in the end.

We evaluate the performance of different DD sequences on various well-known quantum applications, including Bernstein-Vazirani (BV) algorithm, Hidden Shift (HS) algorithm, Quantum Fourier Transform (QFT), and Graph State (GS). The basis quantum circuit structures of these quantum applications are shown in Fig. 2. BV algorithm implements an oracle function f (x), which represents the dot product between x and a secret string s, with the objective of finding s. HS algorithm constructs an oracle that encodes two functions f and g, and there exists a secret string s such that g(x) = f (x + s). The goal of the algorithm is to find s. QFT is the quantum version of discrete Fourier transform, and is the essential part for many other quantum algorithms, such as Shor's algorithm [START_REF] Shor | Polynomial-time algorithms for prime factorization and discrete logarithms on a quantum computer[END_REF], quantum phase estimation algorithm [START_REF] Kitaev | Quantum measurements and the abelian stabilizer problem[END_REF], etc. Graph State is a quantum state prepared based on a graph. Specifically, we can build a Graph State accord-ing to a given hardware topology, where there is an edge on the graph when the two qubits are coupled, and the connection is represented by a CZ gate. It can entangle all the qubits of the device and is important for error correction.

Moreover, we combine different DD sequences with pulse-efficient transpilation technique and apply them to Quantum Approximate Optimzation Algorithm (QAOA) to solve Max-Cut problem. We generate randomly 3regular graphs using ReCirq [START_REF]Quantum AI team and collaborators[END_REF] and random graphs using Networkx with different degrees as our benchmarks. For example, we use 3-reg-4 to represent a 3-regular graph with degree of 4 (qubits) and rand-4-0.5 for a random graph with 4 nodes/qubits and a probability for edge creation of 0.5 (see Fig. 3).

To evaluate the performance of these noise mitigation and pulse-optimization techniques on circuit fidelity, we choose the following metrics for different benchmarks:

• Probability of Successful Trial (PST) [START_REF] Tannu | Not all qubits are created equal: a case for variability-aware policies for nisqera quantum computers[END_REF]. This metric is defined by the ratio of the number of trials that give the expected result to the total number of trials, and higher is better. It is dedicated to the benchmarks with one certain correct result, such as BV algorithm and HS algorithm.

• Jensen-Shannon Divergence (JSD). It is used to measure the similarity between two probability distributions, and lower is better. It is suitable for Graph State circuit and QFT whose output is a distribution.

• Approximation ratio [START_REF] Harrigan | Quantum approximate optimization of non-planar graph problems on a planar superconducting processor[END_REF]. It is specifically designed to evaluate the performance of QAOA circuit and is defined as C /C min , where C is the expectation value obtained by the quantum computer and C min is calculated by the classical NumPyMin-imumEigensolver. We aim at maximizing the approximation ratio and 1 means a perfect solution.

(a) 3-reg-4 (b) rand-4-0.5

FIG. 3: An example of (a) 3-regular graph with degree of 4, (b) random graph with degree of 4 and edge creation probability of 0.5.

IV. EXPERIMENTAL RESULTS

We perform the experiments on various IBM quantum devices with different qubit numbers and quantum volumes (QV) (see Table II). The size of the benchmarks varies depending on the quantum device. When evaluating the performance of DD sequences on different applications including BV, HS, QFT, and GS, the circuit size varies from 2 to 7 for IBM Q 7 Jakarta. Whereas for other devices, the circuit size changes from 2 to 12, since there are too many noise accumulating to obtain meaningful results for benchmarks with more than 12 qubits. For the QAOA experiments aiming at exploring the combination between DD and pulse-efficient optimization technique, the degree of the two types of graphs (3-regular graph and random graph) ranges from 4 to 6 for IBM Q 7 Jakarta, and 4 to 12 for other devices. The QAOA ansatz has one layer with parameters initialized randomly and optimized using COBYLA optimizer on the simulator. We only execute the final ansatz with optimized parameters on the quantum hardware. All the benchmarks are compiled using Qiskit with the highest optimization level and executed 8192 times.

If the DD sequences contain a large number of pulses, such as XY16 and KDD, it might not be possible to insert them for some small benchmarks whose idle time might not be long enough. Moreover, Hahn echo can only be applied if the inverse of the inserted single X or Y gate is able to be absorbed into the neighboring gates to ensure the equivalence of the quantum state. Therefore, we check if the DD pulses are actually inserted to the benchmark for each experiment, and we remove the circuit without any DD pulses inserted from the result.

In order to clearly show the impact of various DD sequences and pulse-efficient optimization method, we use the relative results for the three metrics. The original benchmark with no optimization method applied is marked as the baseline. PST is divided by the baseline and we use JSD results to divide the baseline, so that we can obtain the relative results when applying DD sequences. If the relative result is larger than one, it means that there is an improvement, and larger is better. Whereas for QAOA experiments, we use the difference between approximation ratio with DD sequence or pulse-efficient technique and the baseline as the relative result due to the possible negative value of the approximation ratio. If the difference is larger than zero, it indicates an enhancement and also larger is better. Note that, each experiment has been repeated three times and results with similar trends were obtained.

We first demonstrate applying DD sequences to various quantum applications. The relative PST results for BV, HS circuits, the relative JSD results for QFT and Graph State circuits are shown in Fig. 4, Fig. 5, Fig. 6, and Fig. 7. Hahn echo is not applicable for all the benchmarks and XY16 and KDD are too long to insert for certain small BV circuits.

For BV algorithm (see Fig. 4), all the DD sequences are able to enhance the circuit fidelity, especially XY4 and XY8. The improvement becomes more significant with the increase of benchmark size, and the trend is similar for all the devices. The dramatic drop of fidelity after a certain size might due to the increasing noise. Things become different for HS circuits (see Fig. 5). This type of benchmark seems to be resilient to DD such that there is not much difference after applying each DD sequence, i.e., the maximum relative PST ratio is around 1.1 for all the devices. Whereas for QFT and Graph State circuits (see Fig. 6 and Fig. 7), the implementation of DD is not encouraging. We could barely find a type of DD sequence that always benefits the results. Especially, DD sequences have a higher probability of introducing negative impact for QFT or Graph State circuits with small size.

To conclude, DD can improve the fidelity for certain applications, such as BV algorithm, but not for all the benchmarks. For the DD-beneficial type of application, all the DD sequences are able to improve the fidelity to different extents. Moreover, the DD performance does not vary a lot across different quantum chips and benchmarks.

Second, we show the results of using QAOA to solve Max-Cut problem for different graphs by employing: (1) diverse DD sequences (see Fig. 8), (2) DD sequences + pulse-efficient optimization technique (see Fig. 9). Most of the DD sequences can help enhance the approximate ratio for QAOA, especially CP, CPMG, and UDD X sequences. There is no clear relationship between the level of improvement and the size of the benchmark. Associating pulse-efficient technique with DD sequences can further improve the approximation ratio. The negative impact induced by DD can mostly be canceled due to the pulse-efficient method. For example, the maximum relative approximation ratio for IBM Q 16 Guadalupe grows from 0.246 to 0.37 after applying the pulse-efficient transpilation method, which is an increase of 50%.

One interesting thing is that, even though KDD is a robust DD sequence that is designed to be robust against pulse imperfection and has shown improvements on Rigetti device, its performance is not as good as expected on IBM quantum devices. For all the quantum applications tested including QAOA circuit, there is no growth of circuit fidelity using KDD compared with other DD strategies and actually it decreases the circuit fidelity for most of the cases. It is important to explore different noise mitigation and pulse-level optimization methods, including DD and pulse-efficient optimization technique to better build the quantum circuits. Here, we list several guidelines that can help the community utilize these methods to design circuits with higher fidelity.

• DD cannot always improve the circuit fidelity and it is highly application-dependent. It is recommended to use DD on specific quantum algorithms, including Bernstein-Vazirani algorithm, and QAOA.

• The robust sequence KDD does not work well on IBM quantum devices for most of the quantum algorithms tested.

• Combining DD with pulse-efficient optimization technique is favorable for QAOA circuits and can advance the approximation ratio significantly.

• DD usually has higher probability to have positive impact on circuits with larger depth which might due to a longer idle time.

B. Future work

As the experiments have demonstrated some discouraging effects when inserting DD sequences for certain quantum applications, exploring the hardware physics behind them becomes important. Application-dependent adaptive dynamical decoupling [START_REF] Das | Adapt: Mitigating idling errors in qubits via adaptive dynamical decoupling[END_REF] and variational algorithm-focused VAQEM [START_REF] Ravi | Vaqem: A variational approach to quantum error mitigation[END_REF] have been proposed to provide the most beneficial DD sequences but with large overhead and limited analysis of different DD strategies. A smarter approach for designing appropriate application-oriented DD sequences considering hardware physics is left to future works. Recently, DD has been proven to be able to suppress ZZ-crosstalk for fixed frequency transmon superconducting device [START_REF] Tripathi | Suppression of crosstalk in superconducting qubits using dynamical decoupling[END_REF]. It is interesting to further investigate the performance of DD-based crosstalk suppression at application-level, so that it can contribute to the quantum parallel executions [START_REF] Niu | Enabling multi-programming mechanism for quantum computing in the nisq era[END_REF][START_REF] Niu | How parallel circuit execution can be useful for nisq computing?[END_REF]. Finally, the pulse-level optimization experiments can be delved to other quantum applications to provide more guidelines, such as VQE [START_REF] Kandala | Hardwareefficient variational quantum eigensolver for small molecules and quantum magnets[END_REF], quantum simulation [START_REF] Georgescu | Quantum simulation[END_REF], etc.

VI. CONCLUSION

Today's quantum hardware is prone to noise in the NISQ era. Therefore, circuit optimization and error mitigation approaches are required to increase the output fidelity. In this paper, we focus on two pulse-level circuit optimization methods: dynamical decoupling and pulseefficient optimization technique. First, we implement various DD strategies on IBM quantum devices including non-universal, universal, and robust ones. Second, we apply these DD sequences to several well-known quantum applications, such as QFT and Bernstein-Vazirani algorithm, to evaluate comprehensively the impact of diverse DD techniques on IBM quantum devices with different characteristics. We also merge DD with pulseefficient transpilation method and investigate them on QAOA circuits to solve Max-Cut problem. Based on the experimental results, we found that DD techniques always show positive impact for some benchmarks, while for others, DD demonstrate some discouraging effects. As there is no overhead-free application-oriented DD approach, we provide a list of design guidelines for users to better understand these pulse-level optimization methods and figure out how to improve the circuit design for various quantum applications.

BV3

BV4 BV5 BV6

iθ 2 [ 2 , 0, π 2 n:

 222 cos(θ)σx+sin(θ)σy ] , φ = 0, π repetition time of the basic DD cycle. (θ) φ : a rotation of θ around the axis defined by φ.

  FIG. 1: (a) Original R ZZ gate with θ = π 3 . The circuit is mapped to the first two qubits of IBM Q 27 Toronto. (b) The transpilation result with basis gate sets {R Z , SX, X, CX}, and the duration is 2176 dt. (c) The transpilation result with basis gate sets {R Z , SX, X, R ZX }, and the duration is 1025 dt. We replace R Z ( mπ 2 ) √ XR Z ( nπ 2 ) by U to simplify the circuit.

FIG. 7 :

 7 FIG. 7:Relative JSD results for Graph State circuits on different IBM quantum devices. Higher is better.

FIG. 8 :

 8 FIG.8: Relative approximation ratio of applying DD sequences for QAOA circuit on different IBM quantum devices. Higher is better.

FIG. 9 :

 9 FIG.9: Relative approximation ratio of applying DD sequences + pulse-efficient technique for QAOA circuit on different IBM quantum devices. Higher is better.

TABLE I :

 I A summary of different DD sequences.

	DD sequences	Pulse implementation
	Hahn echo [19]	X or Y
	CP/CPMG [20, 32]	(X) n or (Y ) n
	XY4 [21]	(X -Y -X -Y ) n
	XY8 [33]	(XY 4 -XY 4) n
	XY16 [33]	(XY 8 -XY 8) n
	UDD [34]	(X) n or (Y ) n with space τ where τ = sin 2 ( πj 2n+2 ), j ∈ {1, 2, ...n}
	KDD [22]

TABLE II :

 II The information of different IBM Q devices.

	IBM Q Jakarta Guadalupe Toronto Montreal
	n	7	16	27	27
	QV	16	16	32	128
		V. DISCUSSION
		A. Design guideline

ACKNOWLEDGMENT

This work is funded by the QuantUM Initiative of the Region Occitanie, University of Montpellier and IBM Montpellier. The authors are very grateful to Adrien Suau for the helpful suggestions and feedback on an early version of this manuscript.