Eric Armengaud
email: gmbheric.armengaud@avl.com

Ismar Mustedanagic
email: ismar.mustedanagic@avl.com

Markus Dohr
email: markus.dohr@avl.com

Can Kurtulus
email: can.kurtulus@avl.com

Marco Novaro
email: marco.novaro@ideasandmotion.com

Christoph Gollrad
email: christoph.gollrad@avl.com

Georg Macher
email: georg.macher@avl.com

Avl List

Migration of automotive powertrain control strategies to multi-core computing platforms -lessons learnt on smart BMS

come   L'archive ouverte pluridisciplinaire

Abstract-An important acceptance criteria for electric mobility is the capability to efficiently use the energy stored in the cells of a battery over the vehicle lifetime. The BMS (Battery Management System) plays a central role by estimating the state of charge (current energy available) and state of health (degradation due to ageing effects) of the cells. Improvement of the estimation quality has a direct impact on the battery and thus vehicle range. It is the target of the INCOBAT project to improve the BMS system by means of new electronic components, new control strategies and new development methods in order to achieve cost reduction and performance (driving range) increase. In this context, the introduction of multi-core computing platforms aim at providing more computing resources and additional interfaces to answer the needs of new automotive control strategies with respect to computing performances and connectivity (e.g., connected vehicle, hybrid powertrains). At the same time, the parallel execution, resulting resources and timing conflicts require a paradigm change for the embedded software. Consequently, efficient migration of legacy software on multi-core platform, while guaranteeing at least the same level of integrity and performance as for single cores, is challenging. In this paper, the lessons learnt during the migration of the BMS control strategies to the INCOBAT BMS computing platform will be presented.

Index Terms-Multi Core, Electric Vehicle, Battery Management Systems

I. INTRODUCTION

N recent years, electric mobility has been promoted as the clean and cost-efficient alternative to combustion engines. Although there are already solutions on the market, mass takeup has not yet taken place. There are different challenges that hinder this process from an end user point of view such as costs of the vehicle, driving range, or infrastructure support. Several of these challenges are directly connected to the battery, the central element of the full electric vehicle (FEV). The costs of the battery sum up to 40% of the total costs of a FEV, and the driving range of a FEV is strongly reduced in comparison to the combustion engine.

The aim of INCOBAT 1 (INnovative COst efficient management system for next generation high voltage BATteries, started in October 2013) is to provide innovative and cost efficient battery management systems for next generation HV-batteries. To that end, INCOBAT proposes a platform concept in order to achieve cost reduction, reduced complexity, increased reliability as well as flexibility and higher energy efficiency. Moore's law [START_REF] Moore | Cramming more components onto integrated circuits[END_REF], stating the doubling of the computer capacity every 2 years, is still a strong enabler for this fast function increase and at the same time cost-per function decrease. The current development trend for computing platforms has moved from increasing the frequency of single cores to increasing the parallelism (increasing the number of cores on the same die) to limit the power dissipation while improving the performance. Multi-core and many-core technologies have strong potential to further support the different technology domains, but at the same time present new challenges. Hence, the automotive industry is facing a growing gap between the technologies and required level of expertise to make best use of them. The computing platforms are becoming more and more high-performance with concurrent computing capabilities, larger embedded memories as well as increasing number of integrated peripherals. Low-level mechanisms (e.g., memory protection, diagnostics) typically provided by the basic software or operating system are now being moved into the microcontroller. The complexity of these computing platforms is very high, the related user guides is made of several of thousands of pages. Regarding automotive operating systems and low-level basic software (BSW), the AUTOSAR approach is following a similar trend by standardizing several tens of BSW modules in several tens of thousands pages of specification. Similarly for the application software (ASW, e.g., control strategy for hybrid powertrains), the complexity is already very high and still growing by the introduction of new applications such as advanced driver assistance systems (ADAS) or predictive energy management strategies.

Additionally, the functional integration of the control strategies (e.g., transmission with combustion engine and e-drive) further raises the complexity of the resulting application. The automotive industry is confronted to the central question how to migrate, optimize, and validate a given application (or set of applications) on a given computing platform with a given operating system. A knowledge transfer is required to take over the role of control system integrator and identify the application requirements (both functional and non-functional) and to perform a mapping to the SW and HW architecture. The quality of this mapping has a direct impact on the performance of the control system, and thus of the entire mechatronic system. Main contribution of this paper is to summarize the lessons learnt during the migration of the existing BMS control strategy to the INCOBAT BMS platform based on multi core technology. The paper is organized as follow: Section 2 introduces the INCOBAT project as well as the BMS platform based on Infineon AURIX TM CPU. Section 3 will present the AVL BMS core functions (state of charge SoC, state of health SoH and State of Function SoF estimation) and the adaptation at the functional level that were performed to take advantage of the multi-core platform as well as to enable migration. Section 4 discusses the migration at software integration level and summarizes performance increase achieved. Finally, Section 5 concludes this work.

II. THE INCOBAT PROJECT -AN OVERVIEW

The aim of INCOBAT is to provide innovative and cost efficient battery management systems for next generation HVbatteries. To that end, INCOBAT proposes a platform concept in order to achieve cost reduction, reduced complexity, increased reliability as well as flexibility and higher energy efficiency [START_REF] Armengaud | Improving HV battery efficiency by smart control systems[END_REF]. The targeted outcomes of the project are:  Very tight control of the cell function leading to an increase of the driving range  Radical cost reduction of battery management system  Development of modular concepts for system architecture and partitioning, safety, security, reliability as well as verification and validation, thus enabling efficient integration into different vehicle platforms. To achieve these ambitious targets, the technical approach chosen in INCOBAT primarily relies on the following 12 technical innovations (TI) regrouped into four innovation groups (see Fig. 1):  Customer needs and integration aspects: ensures a correct identification of customer needs and enables an efficient integration into different platforms. This is supported by the use of mission profiles (TI-01)in order to take into account the different driving styles of the customers, the different traffic conditions in the same scenarios and the different tracksand by the integration into a demonstrator vehicle (TI-12)

 Transversal innovation: consistent concept and specification. This second group targets the optimization of the system architecture and its consistent description over the technologies and over the system hierarchies. This aspect aims at providing a consolidated basis in order to simplify later industrialization of the proposed technologies. This includes the TI-02 "Model-based systems engineering" to improve correctness / completeness / consistency of system specification, the TI-03 "System architectureefficient partitioning of the functionalities" for system optimization at BMS or even vehicle level and the TI-04 "Integration of multiple functionalities" to reduce the number of electronic control units (and thus related costs) in the vehicle.  Technology innovation: E/E control system: This third group aims at improving the components of the E/E control system. Regarding the electronic parts, it consists of TI-05 "Multicore computing platform for additional computing resources" and the TI-06 "Smart and integrated module management unit". From the software part, this is achieved by the TI-07 "Modular SW platform" and by TI-08 "Improved BMS control algorithms"  Transversal innovation: improving system maturity: This last group targets the evidences related to the trust on the technical solutions with respect to correct operation (TI-10 "Design and validation plan including reliability consideration"), functional safety and security (TI-09 "Definition and integration of safety and security concept") as well as reliability (TI-11 "Reliability and robustness validation"). This group of technical innovations is an indicator for the maturity of the proposed technology and further provides information on the efforts required for proper integration and validation of the system. The INCOBAT BMS CCU (see Figure 2) is based on the Infineon multicore processor AURIX TC275 with an innovative multicore architecture [START_REF] Brewerton | Practical use of AutoSAR in Safety Critical Automotive Systems[END_REF]. This device supports the concurrent execution of mixed ASIL functions up to ASIL-D [START_REF] Schneider | Multicore vs. Safety[END_REF]. It offers a rich set of peripherals such as A/D converters and timers for data capturing and it has a reasonable number of IOs to support BMS applications. In conjunction with the specific power-supply ASIC TLF35584 it is possible to supply the CCU and support ISO26262 requirements with a minimum number of components. Regarding the SW developments in INCOBAT, a modular development platform is required. Hence, the control strategy and the application software in general are expected to come from different providers and to require different levels of criticality. The activities of SW architectto define the SW blocks as well as their interfacesand the activities of SW integratorintegrating the different SW modules and ensuring correct operation of the entire control systemare especially challenging in the context of automotive supply chain with constraints related to functional safety (ISO 26262 [START_REF]Road vehicles -Functional safety[END_REF]). A modular platform is required to enable the distributed development and flexible deployment of different control strategies and applications in an efficient way. For the SW developments in INCOBAT, the proposed common, modular software development platform consist of:  a layered SW architecture, consisting of several layers and components as well as their interfaces, providing access for the applications to the underlying HW capabilities,  a suitable SW development tool chain, which supports the application SW developers by means of an effective and consistent development process to seamlessly integrate their particular applications to an overall BMS.

III. FUNCTIONAL MIGRATION OF THE BMS CONTROL STRATEGY

A. Model-based battery state estimation 1) Introduction Accurate estimation of battery parameters such as SoC, SoF and SoH requires model-based estimation methods in which a representative model of the battery is utilized as part of the algorithm. Incorporation of a battery model enables many possibilities including:  Capturing expected behavior, to be compared with actual measurements for inference of parameters causing the deviation,  Ability of algorithms to handle different operating conditions and usage scenarios,  Generation of a "prediction error" signal that is necessary for modern estimation methods such as Kalman Filters, or other similar types of Observers.  Ability of the algorithms to be adapted to different cell types, and cell chemistries via only adapting the incorporated battery model  Possibility of making predictions of various battery behaviors as well as battery condition Model based approach to battery state and parameter estimation, therefore has many advantages. An overview of how a battery model may be used as part of the algorithms is given in Figure 3. This architecture uses the output of the battery model for the prediction step of the estimation method, and compares it with the actual measurement to generate an "error" signal. This error signal is then utilized as part of the algorithm to calculate important signals such as SoC and eventually to update the model and adapt it to the actual battery and operation conditions. This architecture is what enables adaptation to various operating conditions and handles situations where the model may not represent battery behavior accurately. The battery model has to be constructed and parameterized in accordance with the type of cell used in the target application since there are large differences of behavior between various types. The model parameterization requires test data that captures the behavior of the cell (e.g. terminal voltage response, surface temperature) under various operation conditions in terms of ambient temperature and usage (i.e. load current). However, there are certain drawbacks included with usage of model based estimation techniques, such as relatively high computation power demand, necessity of a high quality parameterization to obtain satisfactory performance and the necessity to ensure stable behavior over all possible operating conditions.

2)

SoC Estimation The State-of-Charge is defined as the percentage of the maximum possible charge that is present in the battery. The SoC can't be measured directly, but an accurate SoC on pack and on cell level is mandatory for the energy management control system ( The ordinary KF can be used for linear models. Since the battery is a highly nonlinear system, a Kalman filter is necessary, which permits the use of nonlinear models. The Extended Kalman Filter has the ability to handle such kinds of models. A big advantage of a KF is that it considers measurement and process noise due to voltage sensor inaccuracy, temperature fluctuations etc. to estimate the states of the battery. One of the first use of an EKF (Extended Kalman Filter) for SoC estimation of lithium-batteries is published in [START_REF] Plett | Extended Kalman filtering for battery management systems of LiPB-based HEV battery packs Part1. Background[END_REF], [START_REF] Plett | Extended Kalman filtering for battery management systems of LiPB-based HEV battery packs Part 2. Modeling and identification[END_REF] and [START_REF] Plett | Extended Kalman filtering for battery management systems of LiPB-based HEV battery packs Part 3. State and parameter estimation[END_REF]. These papers describe the mathematical background, the modeling of the battery with its identification requirements and the final implementation of the EKF for SoC estimation. Additionally, an algorithm for SoH estimation is presented. The publications [START_REF] Plett | Sigma-point Kalman filtering for battery management systems of LiPB-based HEV battery packs Part 1. Introduction and state estimation[END_REF] and [START_REF] Plett | Sigma-point Kalman filtering for battery management systems of LiPB-based HEV battery packs Part 2: Simultaneous state and parameter estimation[END_REF] describes a modified KF. The Sigma-Point Kalman filter is a more accurate estimation approach, although the computational demand is of the same order as EKF. A comparison between SPKF and EKF show an improvement of the estimation. Moreover, an advanced algorithm is presented, where state and parameter estimation is done simultaneously. A summary of the mentioned algorithms is presented in [START_REF] Plett | Battery management system algorithms for HEV battery state-of-charge and state-of-health estimation[END_REF]. A slightly modified Kalman filter with lower computational demand is used for SoC estimation as part of INCOBAT. The algorithm runs on module level. The SOC-function estimates:  Module States: States of the battery model (model has 6

states, including SoC of the battery)  Module State of Charge  Module OCV: Open circuit voltage of the module The relatively high computational demand makes it very challenging to run the KF on cell-level with state-of-the-art battery management systems (BMS), especially for battery packs intended for high voltage applications. That is why SOC-estimation is done on module level.

3)

SoH Estimation The State-of-Health is a measure of the condition of the battery compared to the fresh battery. It is characterized in the loss in capacity and the increase in resistance. To assure correct SoF and SoH estimation of a battery pack, it is necessary to have information of the cell-SoC and the cellresistance as well. An algorithm to estimate these values is described in [START_REF] Kortschak | Detection Method of Battery Cell Degradation[END_REF]. The basic mathematical background of the estimation approach is described in the next sections. The cell observer is used to determine the deviation of each cell compared to the mean module state already estimated as part of the SoC function and the module resistance observer. The aim is to utilize a much simpler linear algorithm for computational efficiency. The used model is a linearized model that describes the deviation of each cell from the module mean rather than full cell dynamics. The algorithm is a modified version of the Recursive Least Squares, a well-established algorithm commonly used in state and parameter estimation problems.

4) SoF Estimation

The State-of-Function consists of measures for the ability to fulfill the application specific function of the battery. The BMS has to calculate current and voltage limitations such that the battery is operated in a safe operating mode and the performance and lifetime targets are met. State of Function calculation is responsible for determining available functionality of the battery, which would be either current or power that can be supplied to the powertrain, considering the maximum allowable cell voltage and the maximum allowable operating current. The SoF calculation is based on a prediction of future cell voltages for a calibrated prediction time with an electrical model. Basically, the limits are calculated for each cell, as it is shown in Figure 4. The flowchart shows an iteration over the maximum number of cells in the battery pack. For SoF limits, the results of the 'worst' cell (e.g. highest inner resistance) are taken into account. In more detail, the calculation estimates a current limit, based on the average cell voltage in a module. Afterwards, the limits are corrected by the worst cell in a module. The worst cell is identified by the calculated dSoC and dR values.

The algorithm provides the maximum charge and discharge current limits, as well as the maximum charge and discharge power. The power limitations can be calculated by a multiplication of the current and voltage limits.

5) Function inter-dependencies

6)

Functional migration to parallel computing scheme

From the point of view of the model-based battery state estimation, the main focus is on the analysis of the existing algorithms (currently running on single core computing platform) in order to identify possible improvements with respect to SoX estimation accuracy while making use of the additional computing power of a multicore processor. The main target is to improve the accuracy of the estimations from a group of cells (e.g. a module) down to single cell level. In our case, the sensor platform already provides the information required to measure the cells. The challenging factor here is to run a dedicated instance of the existing model estimation algorithm for each single cell instead of one instance for a group of cells. The required computing power (number of algorithm instances running in parallel) is therefore directly dependent on the modelling accuracy (reduction of the number of cells taken into consideration for one algorithm instance). A higher modelling accuracy provide more accurate information on the status of the cells, thus moving the limits for the use of each single cell (and therefore of the entire battery at the end) from a conservative boundary to a more real limit. Consequently, it can be assumed that the range of the vehicle and cycle life of the battery can be increased due to the precise estimation approach.

B. Electrochemical Impedance Spectroscopy (EIS) 1)

Introduction to EIS This methodology measures dielectric properties of a medium as a function of frequency. In particular, as applied to the battery cells, the goal of the EIS is to determine the impedance parameters, and the state of health (SoH) of the cells as a function of the impedance. In order to successfully determine the EIS spectrum it is necessary to take into account certain inherent problems in the method and the component under test. The EIS analysis is based on the following prerequisites: a) the system must be linear b) the system parameters should not vary over time c) the system must be single input, single output (SISO) The lithium battery is not generally satisfying these requirements: therefore, additional assumptions have to be made. First, the characteristic of a battery is not linear: To calculate the impedance it is therefore necessary to proceed to a linearization. The used technique is to identify a working point on the electrical characteristic and to generate a small perturbation of it. Analyzing a small enough portion of a cell's current versus the voltage curve, it is considered to be linear. Therefore, in normal EIS practice, a small (1 to 10 mV) AC signal is applied to the cell: this is small enough to confine the test into a pseudo-linear segment of the cell's current versus voltage curve. Second, the battery parameters are not constant: in general, even with open battery-terminals (i.e. zero current), the battery voltage varies over time depending on the previous history. To allow the stabilization of the battery voltage it is necessary to wait for the conditions of electrochemical balance in the battery. The required time, also referred to as "settling time" or "relaxation time" depends on the temperature (ion mobility) and is estimated in the order of a few hours. A measurement made before reaching the equilibrium condition produces data with variations especially in the lower part of the spectrum. These variations are more or less evident depending on the imbalance inside the cell. Lastly, the voltage in a battery does not depend exclusively on the current flowing through it, but also on other parameters, in particular temperature and SoC. During each EIS measurement, these parameters must remain constant, in order not to influence the output voltage. In general, it must be ensured that the battery open circuit voltage does not vary within the range of the test, or this change will be computed in the spectrum of impedance. Based on the assumptions above, the stimulus signal needed for the EIS test shall have the following characteristics:  The spectrum of the stimulus shall adequately cover the whole frequency range that has to be analyzed (typically from 0.01 Hz to 1 kHz)  The signal amplitude shall be "small enough" to avoid triggering any nonlinear response in the battery As a drawback, the smaller is the signal amplitude, the worse is the signal to noise ratio. For this reason, a dedicated HW solution is needed to obtain a good resolution in the acquired signal. In particular, the proposed solution is:  Amplify and adequately filter the signal, due to the small signal amplitude  Remove the OCV voltage, which is not useful for EIS measurements  Read each voltage value through a differential amplifier All of the above-mentioned features are implemented in a dedicated EIS daughterboard, working together with the CCU-BMS board. In particular, the OCV has to be measured before applying the stimulus signal, and then removed through the dedicated DAC; this is part of the features implemented in the EIS software. Finally, the stimulus signal generated by the daughterboard is a voltage signal; an external amplifier (in particular, a transconductance amplifier) is needed to drive the current that is injected into the cells. The chosen signal, used as a current stimulus, is a sum of sinusoidal current waveforms with predefined frequencies, in the range 0.01 Hz to 1 kHz, with a selectable current amplitude.

2)

The resulting EIS algorithm The EIS algorithm (Figure 7) injects a known current stimulus into the battery cell, reading the resulting voltage. Due to the assumptions described previously, the EIS algorithm will run after a relaxation time, needed for the battery parameters to reach a steady condition. The measured signals from the battery shall then be analyzed, for each frequency, to determine the spectrum of the signal at that frequency. The idea is to correlate the measured signal to the input waveform, to obtain magnitude and phase information about the analyzed signal. The response waveform from the battery typically has a DC offset, harmonic distortion components, and noise components generated by the cell. Nevertheless, the element of the measured signal, which needs to be analyzed, is the one at the same frequency as the generator waveform. All of the spurious components of the measured signal need to be rejected so that accurate measurements of the fundamental signal at the generator frequency can be made. Z(w i )

Figure 7: EIS algorithm overview

The measured system output is multiplied by both the sine and cosine of the test frequency ω. The results of the multiplications are then fed to two identical integrators, where they are averaged over T seconds. As the averaging time increases, the contribution of all unwanted frequency components go to zero and the integrator outputs become constant values which depend only on the gain and phase of the system transfer function at the test frequency.

Harmonics are rejected by the correlation process, and noise is rejected by averaging the signal over a number of cycles; the averaging associated with the correlating frequency response analyzer acts as a band pass filter with center frequency ω. As the average time T increases the bandwidth of the filter becomes narrower, thus the corrupting influence of wide band noise is increasingly filtered out as the correlation time is increased.

Averaging over a complete cycle avoids certain measurement errors associated with offsets on the system output; the performed simulations demonstrates that acquiring on a time window of three complete periods, we obtain an effective rejection of all frequencies above 0.1 Hz. Since the minimum frequency is 0.01 Hz, three complete periods corresponds to 300 seconds. The result of the correlation process is made up of two components one of which is referred to as the Real (or in phase) component, the other is the Imaginary (or quadrature) component. By performing simple mathematical operations on these raw measurement results, it is possible to obtain the magnitude and phase of the impedance.

3) EIS Software implementation

The EIS Software consists of several SW components, using resources either directly from the Aurix microcontroller or through the EIS daughterboard In particular, the EIS consists of:  Complex Device Drivers (CDD): In particular, the application layer implements the EIS algorithm, as described above; the flowchart in Figure 8 shows an overview of the algorithm, that is executed for the whole time of the test (300 seconds), subsequently using then the last output from the integrators to compute the impedance. The resulting SW architecture is depicted in Figure 9.

4)

Functional migration to parallel computing scheme Since the EIS algorithm has been developed from scratch within the INCOBAT project, no specific migration from single core (sequential) to multi-core (parallel) computing scheme was required. The challenges are focused on the proper SW integration and resource managements. This will be discussed further in the next section.

IV. SOFTWARE INTEGRATION ON MULTICORE PLATFORM

A. Software development environment

The large variety of use cases as well as business organization is leading to different requirements on the development framework and build environment:  Flexible configuration of source files, include files and directories for building code for each core. This targets increase in build efficiency as well as constructive integration [START_REF] Kopetz | Real-Time Systems: Design Principles for Distributed Embedded Applications[END_REF] by the capability of updating a core independently from each other  To have a sufficient intellectual property (IP) protection linking of external pre-compiled objects / libraries to the main binary of each core shall be possible. This is especially required in case of distributed developmentmeans different teams / company integrating their IP into a common computing platform  Adaptations to other compilers shall be possible with less effort.  Integration of additional tools shall possible with less effort.

These two last items are related to the distributed development of the entire SW by different teams relying on different development processes and consequently different tools In the context of INCOBAT, the development environment shall be a low cost solution with capabilities to be deployed by each INCOBAT partner while minimizing the licensing costs. To meet all these requirements for the INCOBAT project a set of tools was used to establish the SW development framework. The basic configuration of the SW development environment consists of a standard set of make files and target rules, a common memory mapping and compiler associated make and linker files. To allow utilization of the ERIKA 2 operating system from Evidence into the integration and build process, the command line interface (CLI) from RT Druid was integrated in the build environment in form of make target rules. The OS configuration, including the definition and core allocation of counters, alarms, task, spinlocks and resources is done via OSEK implementation language (OIL) file, which is feed in RT Druid for the generation of the OS Erika related code and header configuration files for each core. The build process is setup in such a manner that for each core a separate binary image is generated. This allows SW updates on one core without the need of rebuilding the other cores. Of course, this mechanism is only applicable if the applied changes do not affect the other's core SW and if the SW of the different cores can access peripherals of the microcontroller only via one dedicated interface. The SW code allocation to the different cores is done statically via one manual configurable make file, in which for each core application source and include files or directories, pre-compiled objects and libraries can be setup. The memory allocation is done aligned to the AUTOSAR memory mapping approach and configured memory sections in the linker script. Depending on the currently identified CPUx in generation, the linker performs allocation of code and data to predefined flash and core local data scratchpad RAM memory sections. The shared data of the cores is defined and allocated by the master core and placed into the local data scratchpad memory of that core which is the producer of the data element. The exchange of the memory information to the slave cores is done by dumping of master core's binary shared memory sections and export to a separate shared sub linker file which is generated during the build of the master core. The sub-linker script is used in the later build phases of slave cores for address resolutions of the shared data elements. Another important aspect for the SW development environment was tool integrationin our case the mapping of the system information with the SW development framework.

During the scope of the project, different tool interfaces for generating AUTOSAR aligned SW information where generated. The proposed tool interfaces mainly relies on four level of exchangesall aligned with the AUTOSAR or OSEK standard. The first level (AUTOSAR tool-bridge) aims at describing the SW components (SW-C) and their interfacesand serves proper integration of control strategies and application SW. The second level (RTE configuration) targets the description of the real-time environment for according configuration. The third level focuses on basic SW (BSW) configuration, while the fourth level aims at describing the operating system (tasks available in the system and their related options). More information is available in [START_REF] Macher | Integration of Heterogeneous Tools to a Seamless Automotive Toolchain[END_REF], [START_REF] Macher | A Model-Based Configuration Approach for Automotive Real-Time Operating Systems[END_REF] and [START_REF] Macher | Automated Generation of Basic Software Configuration of Embedded Systems[END_REF].

B. Software architecture

To ensure modularity and reusability the SW architecture was split in several layers aligned to AUTOSAR:  Infineon iLLDsimilar to AUTOSAR MCAL, providing abstraction to HW I/O's, other peripheral modules, and startup code for the cores  OS Erika -OSEK/VDX certified asymmetric operating system, where each core has its own copy of the OS instance  BSW including complex device drivers and other coded BSW services  ASW and ASWIL, to reduce the complexity of the system each ASW component is accessing its data via separate interfaces from the BSW or IOC module. During the execution of the function group, each core is using its local buffered data wherever possible to minimize execution time caused by inter-core accesses and remote blocking. For the multicore capabilities, several SW functionalities were used similar to the currently defined and supported AUTOSAR concepts:  Synchronized master slave startup and shutdown approach of the cores. During startup of master core 0 the other two slave cores are in idle. They startup with a synchronization barrier during startup of the OS. During shutdown, the reversed order is used and master core waits until synchronized shutdown of slave cores.  Functional based inter-core data exchange of single signals or groups similar to the AUTOSAR Inter-OS-Applicator communicator (IOC).  Usage of spinlocks to guarantee data consistency for coreto-core data exchange. The spinlock mechanism was combined with immediate suspension of interrupts in order to reduce the time of remote blocking. Additionally, to prevent from deadlocks nested acquisition of spinlocks was avoided.

C. Verification environmentmini-HiL

An important target and basis for the proper SW development, verification and integration is the deployment of appropriate test environment. Hence, the test environment shall be flexible While MiL and SiL are state of practice and will not be discussed any further, a dedicated mini-HIL platform for the efficient validation of the BMS (E/E system including satellite units and central controller) has been deployed. This platform is conjointly used for white, grey and black box testingwith different focuses as described previously, For a validation of the SW multicore integration approach and the battery SoX estimation algorithms following HiL test environment was setup, see Figure 12.

In a first phase, the battery estimation algorithms were stimulated with battery cell data via CAN wrapper. For the stimulation of the battery load profiles the test and automation environment from NI Veristand was used in combination CompactRIO HW. With support of the XCP protocol, measurement and calibration access to each core's local memory was established.

D. Performances achieved and lessons learnt

A first comparison of the battery state estimation, while moving the computation accuracy from module to cell level, is shown in the following. The red curves represent the results of the modified estimation approach. Figure 13 illustrates the minimum and maximum cell SoC, while the discharge current limits during the cycle is illustrated in. Figure 14. Especially at the end of the cycle (at 6000s), a higher difference between the estimation on cell-level (red curve) and the estimation on module-level (blue curve) can be recognized. It is important to note that the results achieved in SiL and MiL environment are highly coherent (SoC and SoH difference below 0.05%) therefore confirming the correct integration into the multi-core computing platform. running on each core. Core 0 is managing the BSW and drivers, the battery state estimation algorithm is run on Core 1, and Core 2 is reserved for EIS. It can be noticed that the operating system task is consuming slightly more than 10% core time for this configuration. The computation of the SoX function at pack level requires less than 5% core time. By improving the computation accuracy to cell level, then a computation time of 15% of the core time is required. This illustrates the computing requirement (factor 3) in comparison to computation at module level. At the same time, it illustrates that integration of other functionalities into the AURIX platform (thus reducing the number of electronic control unit) is easily possible. It must be noted that at this stage of the project no particular approach was deployed for the systematic exploration of timing behavior and resource management. The decisions related to startup and shutdown sequences and handling, possible deadlocks, delays caused by remote blocking, memory allocation, as well as ASW functionality scheduling and partitioning were made based on expert knowledge. This step (scheduling analysis) will become essential to ensure an efficient and balanced functionality partitioning and scheduling and is part on ongoing work.

V. CONCLUSION

Accurate understanding and modeling of the physical behaviors of the cells' chemistry is pre-requisite for proper and optimized control of the HV battery, thus moving the limits for use of each single cell (and therefore of the entire battery at the end) from a conservative boundary to a more real limit. Consequently, it can be assumed that the range of the vehicle and cycle life of the battery can be increased due to the precise estimation approach. The continuous advances in chip design (multi-core computing platforms for automotive applications) and embedded SW engineering (AUTOSAR) are providing important basis to deploy complex and accurate control strategies. At the same time, the competences and skills gap is growing apart between the different technologies. The seamless migration of an existing control strategy to a multi-core platform, while considering functional and non-functional requirements (e.g., performances, timeliness, safety), is not an easy task. During this paper, the migration of battery estimation functions (SoX) to an AURIX platform was presented. The migration has led to more accurate battery state estimation and illustrated that the proposed CPU provides enough performances for integration of further functionalities, thus providing the potential for reduction of number of discrete electronic control units within the vehicle. At the same time, an important lesson learnt was the need to proper analyze and manage startup and shutdown sequences, possible deadlocks, delays caused by remote blocking, memory allocation, as well as ASW functionality scheduling and partitioning. These aspects is already part of ongoing work.

Figure 1 :

 1 Figure 1: Technical innovations within INCOBAT

Figure 2 :

 2 Figure 2: INCOBAT BMS CCU Prototype Hardware

Figure 3 :

 3 Figure 3: Overview of Model Based Estimation for Battery States

Figure 4 :

 4 Figure 4: State of Function Flow Chart

Figure 5

 5 provides an overview how the algorithms for SoC, SoH and SoF estimations are implemented. SoC and SoH estimations run on module-level due to the high computational demand, whereas the SOF are calculated per cell. The battery pack SoC is simply the average of the estimated Cell SoC values. To get a plausible SoH of the module resistances of the battery pack, delta SoC and delta R should be considered.

Figure 5 :

 5 Figure 5: Information flow of estimation algorithms

  Figure 6 below shows a functional block diagram of the daughterboard, integrating:  The EIS Command Generator: a voltage DAC needed to generate the stimulus signal  12x EIS cell voltage measurement circuitry (differential amplifier + OCV cancelling + 4th order Bessel antialiasing filter)  2x EIS current measurement circuitry (4th order Bessel antialiasing filter)

Figure 6 :

 6 Figure 6: EIS daughterboard

oFigure 9 :

 9 Figure 8: EIS application flowchart

Figure 10 :

 10 Figure 10: INCOBAT SW development framework

Figure 11 :

 11 Figure 11: INCOBAT SW architecture

Figure 12 :

 12 Figure 12: HiL test environment

Figure 13 :Figure 14 :

 1314 Figure 13: SoC computation at module and cell level

Figure 15 :

 15 Figure 15: Computing resource usage for the three cores

Customer needs and integration aspect 1 Mission profiles Car demonstrator / vehicle validation 12 Transversal innovation: consistent concept and specification bdd [SysML Block Definition] 00_Pow ertrain_Env ironment [BDD_Pow ertrain_Env ironment]

			«block»
			Vehicle
	«block»	«block»	«block»	«block»	«block»	«block»
	Ignition	Steering	Brake	Body	Plug	Vehicle CAN
	«block»		«block»
	Accelerator Pedal	Hybrid Pow ertrain
					«actuator»
					ICE
					«block»
					Electrical Generator
					«block»
					H C U
					«block»
	«actuator»			Transmission
					«block»
					HV Battery

Electrical Traction Motor «block» DCDC Conv erter «block» ACDC Conv erter «block» Pow ertrain CAN ECU «block» LV Battery «block» DCAC Conv erter «block» TCU «block» Gearbox ibd [SysML Internal Block] 00_Pow ertrain_Env ironment [IBD_Pow ertrain_Env ironment] Vehicle

	To Vehicle CAN brake : Brake To Vehicle CAN To Vehicle CAN steering : Steering To Vehicle CAN Ignition Signal ignition : Ignition Ignition Signal From Brake From Steering From Brake From Steering Pedal Position accelerator pedal : Accelerator Pedal Pedal Position Accelerator Pedal Position From Powertrain vehicle bus : Vehicle CAN From Powertrain To Vehicle CAN Ignition Signal hybrid powertrain : Hybrid Powertrain Normal Normal Charging Cable Fast Charging Cable Charging Cable Fast Charging Cable To Vehicle CAN Accelerator Pedal Position Ignition Signal	body : Body Fast Charging Cable Normal Charging Cable plug : Plug Fast Charging Cable Normal Charging Cable	Connection is TBD

  State of Function Calculation). Several methods have been developed in the past:

	 Coulomb counting (Ampere-Hour Counting)
	 Open circuit voltage
	 Neural Networks
	 Heuristic interpretation of measurements, Fuzzy Logic
	 Model-based estimation methods
	o Kalman Filter for battery SoC determination
	o Luenberger Observer
	 Sliding-Mode Observer
	 LPV Observer
	A detailed description of available methods for SoC estimation
	is given in [6]. Furthermore, common algorithms for State-of-
	Health (SoH), State-of-Function (SoF), Remaining useful Life
	(RUL) are introduced.
	Nowadays, the most widely used SoC-estimation algorithms
	are Kalman-filtering techniques. They are based on an
	equivalent circuit model (ECM) of the battery.

http://www.evidence.eu.com/

ACKNOWLEDGMENT

The research leading to these results has received funding from the European Union's Seventh Framework Programme (FP7/2007-2013) under grant agreement n° 608988


J.-F Bercher

C Vignat

A new look at q-exponential distributions via excess statistics

Keywords: q-exponentials, Generalized Pareto Distributions, excess distributions, maximum entropy principle, nonextensive statistics PACS: 02.50.-r, 05.40.-a, 05.90.+m

Q-exponential distributions play an important role in nonextensive statistics. They appear as the canonical distributions, i.e. the maximum generalized q-entropy distributions under mean constraint. Their relevance is also independently justified by their appearance in the theory of superstatistics introduced by Beck and Cohen. In this paper, we provide a third and independent rationale for these distributions. We indicate that q-exponentials are stable by a statistical normalization operation, and that Pickands' extreme values theorem plays the role of a CLT-like theorem in this context. This suggests that q-exponentials can arise in many contexts if the system at hand or the measurement device introduces some threshold. Moreover we give an asymptotic connection between excess distributions and maximum q-entropy. We also highlight the role of Generalized Pareto Distributions in many applications and present several methods for the practical estimation of q-exponentials parameters.

Introduction

Since the pioneering paper by [START_REF] Tsallis | Possible generalization of Boltzmann-Gibbs statistics[END_REF], nonextensive statistics have received an increasing interest in the statistical physics community. As a consequence, the canonical distributions that appear in this context have raised many studies: these distributions are the q-Gaussians under the assumption that the total energy is fixed, and the q-exponentials [START_REF] Tsallis | Possible generalization of Boltzmann-Gibbs statistics[END_REF] when the mean value of the system is fixed. In several cases and a wide variety of fields, experiments, numerical results and analytical derivations fairly agree with the formalism and the description by q-Gaussians or q-exponentials. This includes situations characterized by long-range interactions, long-range memory or space-time multifractal structure. Applications include fully developed turbulence, Levy anomalous diffusion, statistics of cosmic rays, econometry, and many others. Furthermore, recent studies have highlighted their excellent ability to approximate with high accuracy the exact distribution of complex systems [START_REF] Hilhorst | A note on q-gaussians and non-gaussians in statistical mechanics[END_REF]. Let us recall that q-exponential distributions are defined by f q,β (x) = 1 Z q,β e q (-βx)

with e q (-βx) = (1 -β (1 -q) x)

1 1-q , x ∈ [0, A q [ ( 2
)
Email addresses: jf.bercher@esiee.fr (J.-F. Bercher), vignat@univ-mlv.fr (C. Vignat). 1 On sabbatical leave from ESIEE-Paris, France. x e q (-βx) Figure 1. Infinite support q-exponential distributions with β = 1 2-q and q = 1.01 (solid), 1.99 (dash dot) and 1.8 (dash) x e q (-βx) Figure 2. Finite support q-exponential distributions with β = 1 2-q and q = 0 (solid), -∞ (dot), -1 (dash) and 0.7 (dashdot)

where A q is the right endpoint, A q = 1 β(1-q) if q < 1 and A q = +∞ if 1 < q < 2, and where β is a positive scale parameter related to the mean m of f q,β as m = xf q,β (x) dx = 1 β(3-2q) , for q < 3/2. The partition function Z q,β is Z q,β = 1 β (2 -q) , 0 < q < 2.

We will also denote Fq,β (x) the so-called survival function associated to the q-exponential distribution (1):

Fq,β (x) = +∞ x f q,β (z)dz = 1 -F q,β (x) (3)
= (1 -β(1 -q)x) 1 1-q +1 (4)
where F q,β (x) is the probability distribution function.

The q-exponential distributions f q,β (x) are the solutions of the following maximum entropy problem:

max f H q (f ) with +∞ 0 xf q (x) dx = θ and +∞ 0 f (x) dx = 1 (5)
where

H q = 1 1 -q +∞ 0 f q (x) dx -1 (6)
is the Tsallis entropy. We note that the same solution is reached if Tsallis entropy is replaced by using Rényi entropy

S q = 1 1 -q log +∞ 0 f q (x) dx,
since it is a simple monotone transform of the latter. In the limit case q = 1, Tsallis entropy coincides with Shannon entropy and the maximum entropy solution recovers the canonical exponential distribution

f 1,β (x) = β exp (-βx) , x ≥ 0. ( 7
)
The shapes of several q-exponential distributions are shown in Fig. 1 (infinite support case: q > 1 and A q = +∞) and in Fig. 2 (finite support case: q < 1 and A q = 1 β(1-q) ); in all cases, parameter β was chosen as β = 1 2-q so that f q (0) = 1. These figures indicate the high versatility of q-exponentials for modeling probability distributions. An important and alternate characterization of these distributions has been provided by the superstatistics theory as introduced in [START_REF] Beck | Superstatistics[END_REF] (see also [START_REF] Wilk | Interpretation of the nonextensitivity parameter q in some applications of Tsallis statistics and Lévy distributions[END_REF] for related results). In this approach, q-exponential distributions are recognized as scale mixtures of exponential distributions: more precisely, the integral representation

f q, 1 q-1 (x) = +∞ 0 g n (β) f 1,β (x) dβ where g n (β) ∝ β n 2 -1 e -β is a χ 2 n distribution with n = 2
q-1 -2 degrees of freedom, allows to describe a q-exponential random variable as a (q=1) exponential variable with random scale parameter β ∼ χ 2 n . We note that this approach, in a physical setting, recovers similar results by [START_REF] Maguire | The time intervals between industrial accidents[END_REF] [START_REF] Maguire | The time intervals between industrial accidents[END_REF] and [START_REF] Harris | The Pareto distribution as a queue service discipline[END_REF] [START_REF] Harris | The Pareto distribution as a queue service discipline[END_REF]. This superstatistical approach has found many applications in fluctuation theory: for example, in a recent paper, Briggs and Beck [START_REF] Briggs | Modelling train delays with q-exponential functions[END_REF] have highlighted the relevance of the superstatistical approach in the modeling of train delays. Hence, in the one hand, there are two constructions, the maximum q-entropy and the superstatistics one, that exhibit qexponentials. In the other hand, many observations and experimental results assess the relevance of q-exponentials and their ubiquity. We propose in what follows a third and alternate statistical justification of q-exponential distributions, based on the theory of excess values. In section 2, we introduce the "Focus on Excesses" operation, which will reveal as the fundamental operation in our setting. Then, in section 3 we present the Pickands-Balkema-de Haan result as a possible rationale for the ubiquity of q-exponentials in nature. Furthermore, we underline a stability property of q-exponentials, and show that these results form a Central Limit Theorem for the "Focus on Excesses" operation. In section 3.2, we draw attention to the Generalized Pareto Distribution (GPD), which is a q-exponential, and underline its applications outside the statistical physics field. Some examples are worked out in subsection 3.3. In section 4, we address the link between the maximum q-entropy and the distributions of excesses. We draw a possible connection and show that the distribution of excesses converges asymptotically to a maximum q-entropy distribution. Finally, we present and discuss in section 5 the estimation procedures for the parameters of q-exponentials-GPD.

The "Focus on Excesses" operation

Let us consider a one-dimensional positive random variable X (ω) that describes the state of a physical system; we denote as FX its survival function

FX (x) = 1 -F X (x) = Pr {X > x}
with F X (x) the distribution function. Let us also consider a positive parameter u > 0 and assume that we have access only to the values of X that are larger than u : for example, we may have access to the system through a measurement device that, for a bandwidth or quantification reason, does not "see" any value X < u; or we may want to reject all values X < u because we think that they are physically irrelevant. Since parameter u may be large, a wise decision is to shift all these measured values by a factor u : as a consequence, we will in fact record values of the excesses of X over the threshold u, that is values of the conditional random variable X u defined as

X u = X -u|X > u.
We denote X u = F u (X) the transform of X into X u . By Bayes' theorem, the survival function of these excesses writes

FXu (x) = Pr(X > x + u|x > u) = FX (x + u) FX (u) , x ≥ 0. (8)
If X possesses a probability density f X , then X u has itself the density

f Xu (x) = f X (x + u) FX (u) , x ≥ 0. (9)
The transformation of X into X u is illustrated in Fig. 3. The inset depicts the excess distribution, i.e. the normalized tail f Xu (9) of the original distribution f X given in the main part of the figure. For obvious reasons, we call "Focus on Excess" (FoE) operation the transformation of FX into FXu , or equivalently of f X into f Xu . With a slight abuse of notation we also denote F u this transformation. This approach is also known as the "Peaks over Threshold" technique and widely used in fields like climatology and hydrology. A rationale of this method can be found in [START_REF] Leadbetter | On a basis for 'Peaks over Threshold' modeling[END_REF]. The q-exponential distribution f q,β (x) (1) verifies a remarkable stability property with respect to this transformation as described in the following theorem (see [START_REF] Hosking | Parameter and quantile estimation for the Generalized Pareto Distribution[END_REF]). Theorem 1 For any u < A q , the FoE of q-exponential remains q-exponential and leaves the tail index q unchanged:

F u (f q,β (x)) = f q,β ′ (x) with β ′ = β 1-β(1-q)u .
Proof. By equation ( 9), we know that

F u (e q (-βx)) ∝ (1 -β(1 -q)(x + u)) 1 1-q ∝ (1 -β ′ (1 -q)x) 1 1-q , ( 10
) = e q (-β ′ x) for 0 ≤ x ≤ A q -u (11)
where the last term, with

β ′ = β/(1 -β(1 -q)u), is obtained by factoring (1 -β(1 -q)u).
This highlights the particular status of the q-exponential in the FoE transformation. We can also remark that β ′ = β in the special exponential (q = 1) case: the normalized tail of an exponential distribution is exactly invariant by FoE transformation.

The Limit theorem for q-exponentials

Pickands' theorem

An important theorem in probability as well as in statistical physics is the Central Limit theorem; its extension to the nonextensive context has been the subject of several recent papers [START_REF] Umarov | On multivariate generalizations of the q-central limit theorem consistent with nonextensive statistical mechanics[END_REF]. According to this theorem, if a probability density f belongs to the domain of attraction of the Gaussian distribution, and if X i are centered and independently chosen distributed according to f, then the distribution of

Z n = 1 √ n n i=1 X i ( 12
)
converges with n to the Gaussian distribution. Moreover, the Gaussian distribution appears as a fixed point in this context since if all X i are independent and Gaussian in [START_REF] Balkema | Residual life time at great age[END_REF] then Z n is exactly Gaussian for any value of n.

Let us denote by T n the transform T n : f → f n where f n is the distribution of Z n , i.e. the normalized n-fold convolution of f. Consider now a variable submitted to a serie of k successive thresholding operations, with thresholds u 1 , . . . u k : the resulting excess variable, say Z u , is

Z u = F u k • F u k-1 • . . . F u1 (X)
where • is the composition of functions and F ui as defined above. Of course, the excess variable can also be viewed as resulting of a single thresholding operation, with u = k i=1 u i :

Z u = F u (X).
An analogy appears between T n and the FoE transform F u which associates the survival functions:

F u : F → Fu ,
where the real threshold parameter u plays the role of the integer parameter n. Thus, an important question at this point arises: does there exist a limit result, analogous to the CLT, in the context of the FoE transform? Let us first introduce the Fréchet (or heavy-tailed) family of distributions, that will play in extreme statistics the role of the Gaussian domain of attraction in the CLT.

Definition 1 A distribution F belongs the Fréchet domain if it is a heavy tail distribution, that is if its survival distribution F writes F (x) = x -1 γ l (x) , ∀x ≥ 0
where l (x)is a slowly varying function:

lim x→+∞ l (xt) l (x) = 1 ∀t > 0
and parameter γ is called the tail index of F.

Pickands' theorem, which can be viewed as the analogue of the CLT in excess statistics, is the following. Theorem 2 A necessary and sufficient condition for F to belong to the Fréchet domain is

lim u→Aq FXu (x) -Fq,β(u) (x) ∞ = 0
for some function β (u) and with

q = 2γ + 1 γ + 1 . (13)
This result means, roughly speaking, that if F is in the domain of attraction of the Fréchet distribution, i.e. has heavy tails with tail index γ, then its FoE transform converges (in the infinite norm sense) to a limit distribution which is nothing but a q-exponential with q given by (13). To the accumulation of random variates in the classical CLT, with n → ∞, corresponds here an increase in the threshold u in the FoE. In the classical CLT, we know that the limit distribution (the Gaussian) once reached remains stable by convolution of Gaussian pdf (addition of Gaussian variables). Accordingly, theorem 1 shows that q-exponentials are stable by the FoE operation. This means that once the FoE operations have converged to the limit distribution, then it remains stable by further applications of FoE. In other words, q-exponentials are a fixed point of the FoE operation.

Generalized Pareto Distributions

In the original Pickands' formulation [START_REF] Pickands | Statistical inference using extreme order statistics[END_REF] (see also the similar work by Balkema and De Haan [START_REF] Balkema | Residual life time at great age[END_REF]), the name coined for the limit distribution of excesses is Generalized Pareto Distribution (GPD). Its survival function is given by

FX (x) = P (X > x) = 1 + γ σ x -1 γ + for x ≥ 0, (14)
where x + = max (x, 0). Its pdf is

f X (x) = 1 σ 1 + γ σ x -1 γ -1 + for x ≥ 0, (15)
where σ and γ are respectively the scale and shape parameters. A location parameter µ can also be specified according to

f X (x) = 1 σ 1 + γ σ (x -µ) -1 γ -1 + for x ≥ µ.
Obviously, for µ = 0, the GPD coincides with the q-exponential, with

γ = - 1 -q 2 -q and σ = 1 β(2 -q)
The interesting point is that the GPD is employed outside the statistical physics field for modeling heavy tailed distributions. Examples of applications are numerous, ranging from reliability theory [START_REF] Asadi | The mean residual life function of a k-out-of-n structure at the system level[END_REF], traffic in networks [START_REF] Hsieh | An exploratory first step in teletraffic data modeling: evaluation of long-run performance of parameter estimators[END_REF], hydrology [START_REF] Van Montfort | Generalized Pareto Distribution applied to rainfall depths[END_REF][START_REF] Coles | A fully probabilistic approach to extreme rainfall modeling[END_REF], climatology [START_REF] Coles | Extreme value modelling of hurricane wind speeds[END_REF][START_REF] Harris | Generalised Pareto methods for wind extremes. Useful tool or mathematical mirage?[END_REF], geophysics [START_REF] Pisarenko | Characterization of the frequency of extreme events by the Generalized Pareto Distribution[END_REF], materials science [START_REF] Shi | Application of the Generalized Pareto Distribution to the estimation of the size of the maximum inclusion in clean steels[END_REF][START_REF] Fullwood | Generalized Pareto front methods applied to second-order material property closures[END_REF], radar imaging [START_REF] Lacour | Statistical characterization of active sonar reverberation using extreme value theory[END_REF] or actuarial sciences [START_REF] Cebrian | Generalized Pareto fit to the society of actuaries' large claims database[END_REF]. These uses are related to the POT (Peaks over Threshold) method [START_REF] Van Montfort | Generalized Pareto Distribution applied to rainfall depths[END_REF][START_REF] Leadbetter | On a basis for 'Peaks over Threshold' modeling[END_REF][START_REF] Lechner | Assessment of 'Peaks over Threshold' methods for estimating extreme value distribution tails[END_REF], and the underlying rationale is of course the fact that the Balkema-de Haan-Pickands theorem asserts that the distribution of excesses over a threshold often follows approximately a GPD for large values of the threshold.

Some examples

Let us now detail Pickands' theorem for classical distributions that belong to the Fréchet domain, namely: the alphastable, Student and Cauchy distributions.

-The alpha-stable case

In the case of an alpha-stable distribution with parameter 0 < α < 1, we know from [25, Th.1 p. 448] that the tail behavior of the survival distribution is

FZ (z) ∼ 1 Γ (1 -α) z -α so that the tail index is γ = 1 α and FZu (z) ∼ 1 + z u -α
.

We note that the Cauchy distribution (α = 1) can not be deduced from this result.

-The Student-t case

A Student -t distribution with m = 2ν -1 degrees of freedom writes

f X (x) = A ν 1 + x 2 -ν , x ∈ R with A ν = Γ(ν) Γ( 1 2 )Γ(ν-1 2 )
. Thus an equivalent of the Student survival distribution is

FZ (z) = A ν 2ν -1 z -2ν+1 ,
the tail index is γ = (2ν -1) -1 and the excess distribution behaves as

FZu (z) ∼ 1 + z u -2ν+1
.

-The Cauchy case As a special case, the Cauchy distribution corresponds to m = 1 degrees of freedom and

FZu (z) ∼ 1 + z u -1 . (16)
In nonextensive statistics, the Student-t distributions described above are well known under the name of q-Gaussian distributions (q > 1). What is shown here is that the excess distribution of a q-Gaussian distribution with nonextensivity parameter q converges to a q-exponential distribution with nonextensivity parameter

q ′ = 1 + q 2 .
We remark that in limit case where q → 1, the q-Gaussian distribution converges to the classical Gaussian distribution.

In this case we also obtain that q ′ → 1, which means that the excess distribution of a Gaussian converges to the exponential distribution. Indeed, Pickands' theorem extends to the Gumbel family of distributions characterized by an exponential tail: the associated excess pdf is simply the distribution given by (7), which is the limit case of the q-exponential (1) for q → 1.

Numerical Illustration

As a numerical illustration, we present the Cauchy case

f X (x) = 2 π (1 + x 2 )
, x ≥ 0.

The survival function of the excess variable X u is

FXu (x) = FX (x + u) FX (u) = 1 -2 π arctan (x + u) 1 -2 π arctan (u)
.

Using the fact that arctan(t) ≈ π/2 -1/t, for t ≫ 1, we have

FXu (x) ≈ 1 u + x / 1 u = 1 + x u -1
which is ( 16) again. It is the survival function of a q-exponential with entropic index q = 3/2. The excess density can be explicited as

f Xu (x) (x) = 2 π 1 1 -2 π arctan(u) 1 + (x + u) 2
The convergence of the survival distribution to a Generalized Pareto Distribution with γ = 1 as in ( 16) is illustrated in Fig. 4. In this figure, we have plotted the normalized pdfs x → uf Xu (ux) = g u (x), corresponding to the normalized random variable Y u = X u /u, for several values of the threshold u. The figure shows the convergence to the GPD limit g ∞ (x) = 1 (1+x) 2 . We note that Pickands' theorem only ensures convergence of the distribution functions, but the Cauchy case is smooth enough to observe a convergence of the normalized pdf as well. We remark on Fig. 4 that small values of parameter u (typically u = 3) ensure an already accurate approximation of f X∞ . Rate of convergence to the limit distribution are examined in Fig. 5, where we report the quadratic and infinite norms of the approximation error, as a function of u:

g u -g ∞ ∞ and g u -g ∞ 2 2
Note that the rate of convergence of the GPD approximation in Pickands' theorem is rigorously quantified in [START_REF] Worms | Vitesse de convergence de l'approximation de Pareto gà c nà c ralisee de la loi des excès -Rate of convergence for the generalized Pareto approximation of the excesses[END_REF].

A connection with Maximum q-entropy

Since there exists an entropic proof of the classical central limit theorem, it is appealing to look for relationships between Pickands' result and the maximum entropy principle. We show here that a distribution of the excesses over a threshold converges, as this threshold goes to infinity, to a maximum q-entropy solution. Details of similar derivations are given in [START_REF] Bercher | An entropic view of Pickands' theorem[END_REF]. Proposition 3 gives a simplified form of a q-exponential distribution function as the solution of a maximum q-norm problem. Then, we state in theorem 4 that the 1-norm, q-norm, and q-expectation of a suitably normalized version of the excess variable (in the Fréchet domain) converges asymptotically to constant values. Finally, with an appropriate choice of q, this shows that the distribution of excesses is necessarily, asymptotically, the q-exponential solution of the maximum q-norm problem (or equivalently of Tsallis' q-entropy maximization). Proposition 3 With √ 2 > q > 1, consider the class of functions

F = G : R + → R; +∞ 0 zG (z) dz < ∞ .
The function

G * (z) = (1 + z) 1 1-q +1 (17)
is the unique solution of the following maximum norm problem:

max G∈F G q q = max G∈F +∞ 0 G (z) q dz such that +∞ 0 zG (z) q dz = (q -1) 2 (q 2 -q -1)(q 2 -2) and +∞ 0 G (z) dz = 1 -q 2q -3 . (18)
Moreover, the corresponding maximum norm is

G * q q = +∞ 0 G * (z) q dz = 1 -q q 2 -q -1 .
Let us now consider a random variable X that belongs to the Fréchet domain with tail index γ < 1 2 , and a normalized version Y u = X u /g (u) of its excesses, where function g (u) ∼ u, as u → +∞. Then, we have the following results. Theorem 4 With γ < q < 1, as u → +∞, the following asymptotics hold

FYu q q ∼ 1 q γ -1 with asymptotical q-expectation +∞ 0 z FYu (z) q dz ∼ 1 1 -q γ 2 -q γ and asymptotical 1-norm +∞ 0 FYu (z) dz = 1 1 γ -1 .
As a consequence, with the choice of

q = 1 + 2γ 1 + γ or γ = 1 -q q -2
we obtain that the excess distribution of any suitably normalized random variable Y = X/g (u), where X belongs to the Fréchet domain and g (u) ∼ u, has asymptotically same 1-norm, q-norm and q-expectation as the q-exponential distribution [START_REF] Coles | Extreme value modelling of hurricane wind speeds[END_REF]. Since the maximum norm problem [START_REF] Harris | Generalised Pareto methods for wind extremes. Useful tool or mathematical mirage?[END_REF] has a unique solution, we deduce that the maximum q-norm distribution G * and the distribution of excesses F coincide asymptotically as u → +∞. The physical significance behind this result is that it enables to connect the distribution of the excesses over a threshold and a maximum qentropy construction, and thus add a new motivation to the use of Tsallis entropy: q-exponentials, as the Tsallis entropy maximizers, can arise in many contexts in which the system at hand or the measurement device introduces some threshold. Moreover, function G * (17) can now be identified as the survival function associated to the q-exponential density

g * (x) ∝ (1 + x) 1 1-q ,
which is itself the solution of a maximum q-norm problem (5) similar to [START_REF] Harris | Generalised Pareto methods for wind extremes. Useful tool or mathematical mirage?[END_REF]. Hence, the asymptotic density of excesses (if it exists) is a q-exponential density. The addition of a scale parameter β(q -1) simply enables to adjust the mean or variance of the distribution.

Parameters estimation

The statistical analysis of heavy-tailed data requires some caution [START_REF] Clauset | Power-law distributions in empirical data[END_REF]. A recent reference for this topic is [START_REF] Markovich | Nonparametric Analysis of Univariate Heavy-Tailed Data[END_REF].

An important issue for applications of q-exponential-GPD and the assessment of their role in real situations is the estimation of their parameters. As far as the tail exponent α = 1/γ is concerned, several methods exist. The corresponding, and relevant problem, in the nonextensive framework is the estimation of parameter q.

The Hill estimator -The Hill estimator [START_REF] Hill | A simple general approach to inference about the tail of a distribution[END_REF] is given by

γ (k) = 1 k k i=1 log X (n-i+1) -log X (n-k)
where X (1) ≤ X (2) ≤ • • • ≤ X (n) are the order statistics of the sample of size n, and k is a smoothing parameter. This estimator, which can be derived by several rationale (including maximum likelihood) is valid for γ > 0, i.e. q > 1.

Pickands' estimator -Another, more general estimator, is the Pickands' estimator, which is valid for all γ and given by

γ (k) = 1 log 2 log X (n-k+1) -X (n-2k+1) X (n-2k+1
) -X (n-4k+1) Both Hill and Pickands estimators are consistent and asymptotically normal [START_REF] Embrechts | Modelling Extremal Events for Insurance and Finance[END_REF][START_REF] Markovich | Nonparametric Analysis of Univariate Heavy-Tailed Data[END_REF].

Maximum Likelihood -The log-likelihood for a set of n independent and identically distributed observations X i , distributed according to a GPD(γ, σ) given by [START_REF] Van Montfort | Generalized Pareto Distribution applied to rainfall depths[END_REF], is:

log p n (X 1 , X 2 , . . . , X n ) = -n log σ - 1 γ + 1 n i=1 log 1 + γ σ X i .
Then, the Maximum Likelihood procedure gives two equations for the unknown parameters. These relations are implicit but can be solved numerically. We obtain

σ = 1 + γ n n i=1 X i 1 + γ σ X i and γ = 1 n n i=1 log 1 + γ σ X i
More details on Maximum Likelihood estimation of GPD parameters can be found in [START_REF] Grimshaw | Computing maximum likelihood estimates for the Generalized Pareto Distribution[END_REF][START_REF] Johnson | Continuous Univariate Distributions[END_REF]. The Maximum Likelihood for q-exponentials has also be discussed recently in [START_REF] Shalizi | Maximum likelihood estimation for q-exponential (Tsallis) distributions[END_REF]. It shall be mentioned that the ML does not exist for γ < -1, which corresponds to the case q > 2, where e q (x) is not normalizable.

Method of moments -When both the mean m X and variance s 2 X exist, that is for γ < 0.5 (q > 0.5), one can apply a method of moments (MOM) since there is a simple relationship between these moments and the parameters γ and σ, namely σ = 0.5m X 1 + (m X /s X ) 2 (19) and γ = 0.5 1 -(m X /s X ) 2 (20) Then it suffices to use some estimates of the mean and standard deviation in [START_REF] Pisarenko | Characterization of the frequency of extreme events by the Generalized Pareto Distribution[END_REF][START_REF] Shi | Application of the Generalized Pareto Distribution to the estimation of the size of the maximum inclusion in clean steels[END_REF] to obtain a (rough) estimate of the GPD parameters.

Conditional mean exceedance -Variations on the method of moments include the so-called "conditional mean exceedance method" which relies on the analysis of the plot of the mean of excesses over u, which is linear according to

E[X -u|X > u] = σ + γu 1 -γ ,
with slope γ/(1 -γ) and intercept σ/(1 -γ). For several values of the threshold u, and therefore estimates of the mean of excesses, a least-squares procedure can be used in order to identify the parameters of the linear model.

Probability Weighted Moments -An alternative method to the MOM, known as "Probability Weighted Moments" (PWM) was introduced by [START_REF] Hosking | Parameter and quantile estimation for the Generalized Pareto Distribution[END_REF][START_REF] Hosking | Some theory and practical uses of trimmed l-moments[END_REF]. The paper [START_REF] Hosking | Parameter and quantile estimation for the Generalized Pareto Distribution[END_REF] also provides and interesting comparison of ML, MOM and PWM estimates. This procedure relies on the definition of a "weighted moment":

m p = E X FX (X) p
where FX is the survival function of the GPD given in [START_REF] Hsieh | An exploratory first step in teletraffic data modeling: evaluation of long-run performance of parameter estimators[END_REF]. For γ < 1 (q > 0), we readily obtain that

m p = σ (p + 1)(p + 1 -γ) .
Then, two weighted moments are sufficient to exhibit the values of the parameters γ and σ; e.g. with m 1 and m 0 , one has γ = 4m 1 -m 0 2m 1 -m 0 and β = 2m 1 m 0 m 0 -2m 1 With more moments, estimated classically, one can look for a least-squares solution. PWM is reported to compete with Maximum Likelihood.

Synthesis and Software -Another method is the Elemental Percentile Method, discussed in [START_REF] Castillo | Fitting the Generalized Pareto Distribution to data[END_REF]. A synthesis of estimation methods of GPD parameters can be found in [33, chapter 20, pages 614-620] and in [START_REF] Beirlant | Statistics of Extremes: Theory and Applications[END_REF].

It shall be mentioned that statistical softwares include estimation methods for GPD-q-exponentials parameters. An example is the "Vector Generalized Additive Models" package [START_REF] Yee | Vector Generalized Additive Models[END_REF][START_REF] Yee | The VGAM package[END_REF]. Another example is the statistical toolbox of the highly employed Matlab TM software, which includes fitting tools for the GPD since its 5.1 version (2006).

Conclusion

In this paper, we have proposed a possible rationale for understanding the ubiquity of q-exponential distributions in nature. The point is that, given an heavy-tailed distributed system, q-exponentials can occur as soon as the measurement device, or the system at hand, involves some threshold. This is the essence of Pickands' result. Furthermore, the q-exponentials are stable by thresholding. This shows that q-exponentials are the limit distributions in a Central Limit like theorem, where the underlying operation is the "Focus on Excesses" operation we have introduced. We have also underlined that q-exponentials can also be recognized as Generalized Pareto Distributions which are of importance in many applications outside the statistical physics field. We have drawn a connection between excess distributions and the maximum q-entropy, by showing that a distribution of excesses converges to a maximum q-entropy solution. Finally, we have presented some procedures for the crucial problem of estimation of q-exponential parameters. We have also provided a large number of bibliographic entries on these different topics. In our view, the rationale proposed here is complementary to the classical maximum q-entropy approach, and to the Superstatistics approach. We hope that it will prompt some new viewpoints in physical applications.

Figure 3 .

 3 Figure 3. Illustration of the FoE transformation

Figure 4 .Figure 5 .

 45 Figure 4. Illustration of Pickands' theorem in the Cauchy case. The figure presents the normalized excess pdfs gu(x), for several values of the threshold u. It shows that the normalized pdf quickly converges, as u increases, to the limit distribution 1/(1 + x) 2 .


Jean-François Bercher
email: bercherj@esiee.fr

Christophe Vignat
email: vignat@univ-mlv.fr

An Entropic View of Pickands' Theorem

It is shown that distributions arising in Rényi-Tsallis maximum entropy setting are related to the Generalized Pareto Distributions (GPD) that are widely used for modeling the tails of distributions. The relevance of such modelization, as well as the ubiquity of GPD in practical situations follows from Balkema-De Haan-Pickands theorem on the distribution of excesses (over a high threshold). We provide an entropic view of this result, by showing that the distribution of a suitably normalized excess variable converges to the solution of a maximum Tsallis entropy, which is the GPD. This highlights the relevance of the so-called 'Tsallis' distributions in many applications as well as some relevance to the use of the corresponding entropy.

I. INTRODUCTION

Generalized Pareto Distributions (GPD) are widely used in practice for modeling the tails of distributions. The underlying rationale is the Balkema-De Haan-Pickands theorem [?], [?], which asserts that the distribution function of the excess variable X -u|X > u (that is the distribution of the (shifted) variable X exceeding a threshold u) converges to a GPD with survival function, as u ∞:

S X (x) = P r (X > x) = 1 + γ σ x -1 γ , (1)
where σ is a scale parameter and γ a shape parameter. For γ = 0, the GPD reduces to the exponential distribution. In applied fields, GPD have encountered a large success since they were obtained as the maximizers of a special entropy, the Tsallis (Harva-Charvat-Daróvczy) entropy [?], with suitable constraints. It is worth mentioning that monotonous transforms of the latter entropy also exhibit the same GPD maximizers: an important example is Rényi entropy [?]. Of course, in a wide variety of fields, experiments, numerical results and analytical derivations fairly agree with the description by distribution (1). This distribution is of very high interest in many physical systems, since it enables to model powerlaw phenomena. Indeed, power-laws are especially interesting since they appear widely in physics, biology, economy, and many other fields [?].

In this communication, we give an interpretation of Pickands' result which links it to the maximum (Rényi/Tsallis) entropy setting with emphasis on the the importance of GPD for modeling excess distributions. This view gives a possible interpretation for the ubiquity of 'Tsallis' (GPD) distributions in physics applications, as well as in other fields, an an argument in support to the use of Rényi/Tsallis entropies.

In the following, we deal with univariate distributions defined on R or on a subset of R. Our approach is as follows. First, we show that the GPD can be obtained as the solution of as resulting of a Rényi-Tsallis maximum entropy problem, with normalization and moment constraints. Second, we consider distributions in the domain of attraction of Fréchet distributions. This family includes for instance Cauchy, Student and Pareto distributions. We characterize the associated q-norm and first moment of the survival function associated to the excess variable X u|X > u. Using an appropriate normalization, we define a variable whose survival function qnorm and moment converge to constant values. Third, we show that it is possible to exhibit a parameter q, as a function of the the parameter of the Fréchet family, such that the maximum entropy and values of the constraints coincide with those of the normalized variable from the Fréchet family. Therefore, since the maximum entropy with given constraints is unique, we obtain that the excess variable from a distribution in the domain of attraction of Fréchet distribution asymptotically follows a Generalized Pareto Distribution (1).

II. SOLUTION TO THE MAXIMIZATION OF TSALLIS' ENTROPY

We first begin by the expression of the solution to the maximization of Tsallis' entropy subject to normalization and moment constraints.

Theorem 1: Consider the set F = {G : R + → R} The maximum Tsallis entropy problem (or equivalently the maximum q-norm problem), with q < 1, defined by

max G∈F S q (G) = max G∈F 1 1 -q   D G(z) q dz -1   subject to +∞ 0 zG(z)dz = µ and +∞ 0 G(z)dz = θ (2)
has for unique solution

G * (z) = α 1 q-1 1 + β α z 1 q-1 (3) with, for q ≥ 1/2 µ = (q -1) 2 q (2q -1) α 2q-1 q-1 β 2 , θ = α q q-1 β (1 -q) q (4)
and

S q (G * ) = α 2q-1 q-1 β (1 -q) (2q -1) . ( 5
)
The mean is not defined for q < 1/2. Proof: We follow here the approach of [?]. Consider the functional Bregman divergence:

B(f, g) = d(f, g)dx (6) = -f (x) q -g(x) q -q (f (x) -g(x)) g(x) q-1 dx (7)
associated to the (pointwise) Bregman divergence d(f, g) built upon the strictly convex function -x q for q ∈ (0, 1). Then let us evaluate the divergence between the distribution G * (z) in ( 3) and any distribution G(z), with G dominated by G * , G(z) ≪ G * (z), and satisfying (2):

B(G, G * ) = - S G(z) q -G * (z) q -α(G(z)G * (z) q-1 -G * (z) q )dz (8) = - S G(z) q dz + S G * (z) q dz, (9)
where S denotes the support of G * (z). The second line follows from the fact that since G and G * both satisfy (2), then, using (3) it is easy to check that

S G(x)G * (x) q-1 dx = S G * (x) q dx.
The Bregman divergence B(G, G * ) being always positive and equal to zero if and only if G = G * , the equality (9) implies that, for q ∈ (0, 1),

S q (G * ) ≥ S q (G) (10)
which means that G * is the distribution with maximum Rényi-Tsallis entropy, with q ∈ (0, 1), in the set of all distributions G ≪ G * satisfying the constraints (2). Values of the constraints (4) and of the maximum entropy (5) follow by direct calculation.

III. THE DISTRIBUTION OF EXCESSES FOR DISTRIBUTIONS

IN FRÉCHET FAMILY In the following, we consider the set of distributions which belongs to the Fréchet domain of attraction. This is the set of all distributions F such that if variables X i are independent and identically distributed according to a distribution in F, then max i=1..N X i converges to the Fréchet distribution. This family typically represents heavy-tailed distributions whose tail behave as a power-law. It was shown by Gnedenko [?] that a necessary and sufficient condition for a distribution to be in the Fréchet domain of attraction is that its survival function satisfies lim z→+∞

S(z) S(cz)

= c a , with c > 0, a > 0. Equivalently, this can also read

lim z→+∞ S(z) = z -a l(z),
where l(z) is a slowly varying function, i.e. a function such that lim z→+∞ l(zt)

l(z) = 1, ∀t > 0. Let us consider the excess variable X u = X -u|X > u. Its survival function is S Xu (z) = S X (z + u) S X (u) .
Theorem 2: Suppose that X u belongs to the Fréchet domain, with S Xu (z) ∼ z -a l (z) , then its q-norm is asymptotically

S Xu q ∼ u aq -1
and its first moment is asymptotically, with a ≥ 2,

+∞ 0 zS Xu (z) dz = u 2 (1 -a) (2 -a)
.

Proof: the q-norm writes

S Xu q = +∞ 0 S X (z + u) S X (u) q dz = +∞ u S X (z) S X (u) q dz = u +∞ 1 S X (wu) S X (u) q dw ∼ u +∞ 1 (uw) -aq u -a l (wu) l (u) q dw ∼ u +∞ 1 (uw)
-aq

u -a dw = u +∞ 1 w -aq dw = u aq -1 , with 1 -aq ≤ 0, since a ≥ 2, q ≥ 1/2.
Of course, we immediately obtain, taking q = 1, that

S Xu 1 = u a -1 .
Similarly, the first moment is

+∞ 0 zS Xu (z) dz = +∞ 0 z S X (z + u) S X (u) dz = +∞ u (z -u) S X (z) S X (u) dz = +∞ 1 u (w -1) S X (wu) S X (u) udw ∼ u 2 +∞ 1 (w -1) w -a dw = u 2 (1 -a) (2 -a)
.

We have a simple corollary to this theorem:

Corollary 1: The survival function S Yu of random variable Y = X/g (u)
, where function g is such that g (u) ∼ u, has asymptotical norms

S Yu q ∼ 1 aq -1 and S Yu 1 = 1 a -1 .
and an asymptotical first moment

+∞ 0 zS Yu (z) dz ∼ 1 (1 -a) (2 -a)
.

Proof: The results for S Yu follow directly from Theorem 2, with S Yu (z) = S Xu (zg (u)) .

IV. THE ENTROPY OF THE DISTRIBUTION OF EXCESSES

Of course, coincidence of the survival function S Yu (z) of the normalized excess Y u with the maximum entropy solution (3) imposes a = 1 1 -q .

Then, it is easy to check that Theorem 3: if X belongs to the Fréchet family, then choosing q < 1 such that a = 1 1 -q , the excess distribution of Y u reaches asymptotically the maximum q-norm solution under constraints asymptotically equal to µ and θ provided α = β = 1.

Proof:

Choosing a = 1 1-q yields S Yu q ∼ 1 aq -1 = 1 -q 2q -1 , S Yu 1 ∼ 1 a -1 = 1 -q q and +∞ 0 zS Yu (z) dz ∼ 1 (1 -a) (2 -a) = (q -1) 2 q (2q -1)
which coincide with the unique maximum q-norm function with constraints µ and θ if and only if α = β = 1.

We can now return to the distribution of the excess without normalization. Then we obtain a Generalized Pareto Distribution, whose shape parameter is given by the exponent parameter of the Fréchet family, and where the scale parameter is a function of the threshold u and of the Fréchet parameter a. This can be stated as follows.

Theorem 4: If X belongs to the Fréchet family, then the survival function of its excesses over a threshold u converges (in the q-norm sense) to

G (z) = 1 + γ σ x -1 γ (12)
with γ = 1/a and σ = g(u)/a for any positive function g such that g (u) ∼ u.

Proof: The proof is immediate by (11). At this point it is still important to emphasize that the GPD (12) enjoys a threshold stability property: the distribution of the excesses over a threshold of GPD remains a GPD, with the same exponent but a different shape parameter. This stability result is certainly also a reason of the ubiquity of GPD in many applications.

Theorem 5: a GPD with parameters γ, σ the distribution of excesses remains a GPD with parameters γ and σ ′ = σ 1 + γ σ u . Proof: As usual, let u denotes the threshold, S X the survival function of the original GPD and S Xu the survival function of variable X u .Then, S Xu (x) = S X (x + u) S X (u) = 1 + γ σ (x + u)

-1 γ 1 + γ σ u -1 γ = 1 + γ σ ′ x -1 γ
where the last expression is obtained after having factored the term 1 + γ σ u -1 γ in the numerator.

We have shown that the distribution of excesses in the Fréchet family, which is a Generalized Pareto Distribution, can be interpreted as a maximum q-entropy (or maximum qnorm) solution. With this result, it is possible to connect the ubiquity of heavy-tailed distributions in physics, economics or signal processing, the distribution of the excesses over a (often unknown) threshold, and a maximum entropy construction. The final paper will present a numerical illustration in the Student case, and will include a simple extension to the Weibull family, wich leads to distributions with exponential tails.


J.-F Bercher

Maximum entropy with fluctuating constraints The example of K-distributions

Keywords: maximum entropy principle, superstatistics, K-distributions PACS: 05.40.-a, 02.50.-r, 05.90.+m

Maximum entropy with fluctuating constraints: The example of K-distributions

Introduction

The technique of compounding, or mixing, is a standard technique for the interpretation or construction of distributions [1, chapter 8]. A random variable may result of the mixing of different populations, or the source of the variable may be unknown and the possible sources characterized by a probability density function. This technique has been employed in [START_REF] Wilk | Interpretation of the nonextensitivity parameter q in some applications of Tsallis statistics and Lévy distributions[END_REF][START_REF] Beck | Superstatistics[END_REF] for recovering the Tsallis (Generalized Pareto) distribution and generalized as 'Superstatistics'. This concept that has met interest [START_REF] Sattin | Multiparameter generalization of nonextensive statistical mechanics[END_REF][START_REF] Bashkirov | On maximum entropy principle, superstatistics, power-law distribution and Rényi parameter[END_REF][START_REF] Queirós | On the emergence of a generalised Gamma distribution. Application to traded volume in financial markets[END_REF][START_REF] Yamano | Thermodynamical and informational structure of superstatistics[END_REF][START_REF] Wilk | Fluctuations, correlations and the nonextensivity[END_REF][START_REF] Abe | Superstatistics, thermodynamics, and fluctuations[END_REF] consists in varying the natural parameter β of a Boltzmann distribution according to some density f β (β) for β. In particular, the Tsallis distribution is associated with a Gamma distribution. Interestingly, this distribution had already been derived in such a way in [START_REF] Maguire | The time intervals between industrial accidents[END_REF]. The objective of this Letter is twofold. First, we relate fluctuations of the intensive parameter β to fluctuations of an observable Ē, and indicate that the relevant distribution of the underlying system can be viewed as the maximum entropy solution with randomized constraints. Second, we work out the case of a Gamma distribution for the observable and draw attention to the physically relevant solution, a K-distribution, that is characterized and discussed.

Email address: jf.bercher@esiee.fr (J.-F. Bercher). 1 On sabbatical leave from ESIEE-Paris, France.

Preprint submitted to Physics Letters A Typesetted 15 May 2008

Consider a Gibbs ensemble where each system of the ensemble is composed of equilibrated subsystems, cells, each of them being described by the standard canonical factor e -βE . The corresponding distribution can be regarded as the distribution that maximizes the Shannon-Boltzmann entropy

S = -f E (x) log f E (x)dΓ (1)
subject to normalization and to the observation constraint Ē = xf E (x)dΓ, and where dΓ is the infinitesimal volume in the phase space. Resolution leads to f E (x| Ē) = e -βx-log Z(β) , where Z(β) = e -βE dΓ is the partition function and where β, the Lagrange multiplier associated to the observation constraint is solution of

S( Ē) = max β {-β Ē -log Z(β)}. (2)
This relation simply indicates that the entropy S( Ē) and the potential Φ(β) = log Z(β) are conjugated functionals while Ē and β are conjugated variables. This also gives the standard relationship

Ē = - dΦ(β) dβ = - d log Z(β) dβ . (3)
Hence, it is stressed that these variables are not independent, and that variations of the thermodynamic β necessarily implies variations of its dual variable Ē and reciprocally. Thus, if the mean energy varies from subsystem to subsystem (from cell to cell), and if these variations are modeled using a statistical distribution f Ē ( Ē), then the global distribution of the system under consideration will be

P (E) = +∞ 0 f E (E| Ē)f Ē ( Ē)d Ē, (4)
where f E (x| Ē) is the Boltzmann distribution with mean Ē. Given a model f Ē ( Ē) of fluctuations of Ē and the relationship (3) linking β and Ē, we readily obtain the distribution for β:

f β (β) = f Ē -dΦ(β) dβ d 2 Φ(β) dβ 2 . (5)
Then, the system distribution can also be written

P (E) = +∞ 0 f E (E|β)f β (β)dβ. (6)
Of course this is similar to the 'Superstatistics' of Beck and Cohen [START_REF] Beck | Superstatistics[END_REF], but with the supplementary ingredient of the relationship between the intensive parameter β and its dual variable Ē, so that a model of fluctuations can be naturally introduced for Ē. At this point, it can also be noted that the dimension of β is the inverse of an energy E.

Illustration in the canonical case

Let us now illustrate this point of view with a very simple example: the probability distributions of velocity and kinetic energy for a macroscopic body consisting of a large number of atoms. In the Gibbs distribution, the energy E(p, q) can be decomposed in the sum of the kinetic and potential energy E(p, q) = K(p) + U (q), so that e -βE(p,q) dpdq = e -βK(p) e -βU(q) dpdq, with p the momenta and q the coordinates. Then the probability for the momenta reduces to Ae -βK(p) dp. We know that for an atom of mass m the kinetic energy is K(p) = p 2 x + p 2 y + p 2 z /2m, where p x , p y , p z are the coordinates of the momentum, and one obtains the classical Maxwellian distribution. Expressing the distribution in terms of the velocities, with p = mv and

v 2 = v 2 x + v 2 y + v 2 z , we have f V (v|β) = 4π βm 2π 3/2 v 2 e -βmv 2 /2 (7) f E (E|β) = 2 √ π β 3 2 √ Ee -βE (8)
for the magnitude of the kinetic energy, which is a χ 2 distribution with three degrees of freedom. If the motion was restricted to a single axis, we would have a single degree of freedom and an exponent 1/2 for β instead of 3/2 in the previous relations. The partition function is given by

Z(β) = +∞ 0 e -βK(p) dp = π 2βm 3 2 , (9)
and we immediately obtain

Ē = - d dβ log Z(β) = 3 2β . ( 10
)
Therefore, in the canonical case, the mean value of the kinetic energy is simply proportional to the inverse of the parameter β. Fluctuations of this mean value Ē can be modeled through a distribution f Ē .

The example of K-distributions

It is well known, e.g. [START_REF] Wilk | Interpretation of the nonextensitivity parameter q in some applications of Tsallis statistics and Lévy distributions[END_REF][START_REF] Maguire | The time intervals between industrial accidents[END_REF], that Gamma fluctuations of the 'inverse temperature' β leads to a Tsallis distribution. Because of the relationship [START_REF] Maguire | The time intervals between industrial accidents[END_REF], this corresponds to an inverse-Gamma distribution for Ē. Another reasonable model for fluctuations of the mean energy Ē is a Gamma distribution since it proves highly flexible for modeling the distribution of positive variables. It also appears as the distribution of the sum of squared Gaussian variables.

Derivation

Therefore, according to [START_REF] Maguire | The time intervals between industrial accidents[END_REF], β is distributed according to the inverse-Gamma distribution:

f β (β) = b a Γ(a) β -a-1 e -b β
. Then, P (E) is given by

P (E) = 2 √ π b a Γ(a) +∞ 0 β 3 2 √ Ee -βE β -a-1 e -b β dβ. (11)
Formula [3.471.9] in [START_REF] Gradshteyn | Table of Integrals, Series and Products[END_REF] enables to perform the integration with respect to β. This formula reads +∞ 0

x ν-1 e -β x -γx dx = 2 β γ ν 2 K ν (2 √ βγ)
with β, γ > 0, and where K ν is a modified Bessel function the of second kind of order ν. This leads to

P (E) = 2 E 1 Γ (a) Γ (3/2) (bE) 3/2+a 2 K -3 2 +a 2 √ bE . ( 12
)
As far as the velocity is concerned, the variable change E = 1 2 mv 2 gives

P V (v) = 4 Γ(3/2)Γ(a) v a+ 1 2 c a+ 3 2 K a-3 2 (2cv) (13)
with c = bm/2.

Main properties and characterizations

These distributions are known as K-distributions. They are extensively used for modeling the statistics of interferences of radiowaves, radar clutter, optical scintillation, ultrasound scattering, etc [START_REF] Jakeman | Significance of K distributions in scattering experiments[END_REF][START_REF] Jakeman | On the statistics of K-distributed noise[END_REF][START_REF] Raghavan | A method for estimating parameters of K-distributed clutter[END_REF][START_REF] Abdi | K distribution: an appropriate substitute for Rayleigh-lognormal distribution in fading-shadowing wireless channels[END_REF][START_REF] Iskander | A method for estimating the parameters of the K-distribution[END_REF][START_REF] Lacour | Statistical characterization of active sonar reverberation using extreme value theory[END_REF][START_REF] Ward | Sea Clutter: Scattering, the K Distribution And Radar Performance[END_REF][START_REF] Weng | Ultrasound speckle analysis based on the K-distribution[END_REF]. Interestingly, a recent application of inverse-chi-square superstatistics in medical statistics has been proposed [START_REF] Chen | A superstatistical model of metastasis and cancer survival[END_REF] with excellent agreement with real data. The shape of distribution [START_REF] Jakeman | Significance of K distributions in scattering experiments[END_REF] and its evolution when the parameter a varies is given in Fig. 1. This distribution is derived in radar imaging as the mixture of a Rayleigh and Gamma distribution [START_REF] Ward | Sea Clutter: Scattering, the K Distribution And Radar Performance[END_REF] and competes with log-normal e.g. [START_REF] Abdi | K distribution: an appropriate substitute for Rayleigh-lognormal distribution in fading-shadowing wireless channels[END_REF] or generalized Pareto (Tsallis) distributions e.g. [START_REF] Lacour | Statistical characterization of active sonar reverberation using extreme value theory[END_REF] for fitting observed data. It is worth mentioning that the estimation of the parameters of this distribution is well documented, e.g. [START_REF] Raghavan | A method for estimating parameters of K-distributed clutter[END_REF][START_REF] Iskander | A method for estimating the parameters of the K-distribution[END_REF][START_REF] Ward | Sea Clutter: Scattering, the K Distribution And Radar Performance[END_REF]. 12) with a = 0.1 to a = 10, and with b = 3a/2 chosen so that the expectation of all distributions is equal to 1.

The case of general mixture of Gamma distributions was considered in [START_REF] Bhattacharya | A modified Bessel function model in life testing[END_REF] where is given a formula similar to [START_REF] Jakeman | Significance of K distributions in scattering experiments[END_REF]. A K-distribution was also obtained in [START_REF] Sattin | Multiparameter generalization of nonextensive statistical mechanics[END_REF] in the limit case of a ratio of two Gamma distributions (that is a F-distribution) for describing fluctuations of β. A model of inverse chi-square for β was mentioned and its asymptotics predicted in [START_REF] Touchette | Asymptotics of superstatistics[END_REF]. These asymptotics are easy to observe: formula [9.7.2] in [START_REF] Abramowitz | Handbook of Mathematical Functions: with Formulas, Graphs, and Mathematical Tables[END_REF] gives the limit form for the Bessel K function for large arguments: K ν (z) ≈ π 2z e -z . Consequently, P (E) behaves as

P (E) ∼ 2 Γ (a) (bE) a-1 2 e -2 √ bE , (14)
and presents an exponential tail. The survival function F (x) = Pr(X ≥ x) then has the limit form F (x) = For small values of the argument, [START_REF] Abramowitz | Handbook of Mathematical Functions: with Formulas, Graphs, and Mathematical Tables[END_REF] again gives a limiting form, formula [9.6.9]: K ν (z) ≈ 1 2 Γ(ν) 1 2 z -ν for ν > 0. For a > 3/2 we then have

P (E) ∼ Γ(a -3 2 ) Γ (a) Γ (3/2) b 3 2 √ E. (15)
The case a < 3/2 follows using the fact that K -ν (z) = K ν (z):

P (E) ∼ Γ( 3 2 -a) Γ (a) Γ (3/2) b a E a-1 . (16)
This last expression shows that the K-distribution presents an integrable singularity at the origin for a < 1, and explains the behavior shown in Fig. 1. In such a case, there is a preferential weighting of small values interrupted by some rare large events, an interesting feature encountered in many experiments. A typical realization obtained with a = 0.5 is reported in Fig. 2. As a final comment, let us recall that the distribution of the ratio of a Gaussian variable and a Gamma variable is a Tsallis distribution. Along the same line, a stochastic model of K-distributed variates can be a simple multiplicative model. Indeed, if z is the product of two Gamma variables x and y, the distribution of z is p Z (z) = p X (x)p Y (z/x)dx/x, which reduces to a K-distribution.

Conclusion

This Letter presents two contributions. First, we have highlighted the fact that within a maximum entropy setting, the β parameter and observable are conjugated, and linked by relation [START_REF] Beck | Superstatistics[END_REF]. Therefore, their fluctuations are also linked. This gives an alternate viewpoint to 'Superstatistics', since the model of fluctuations can be introduced for the observable rather than for the intensive parameter. As an illustration and second contribution, we worked out the natural example of a Gamma distribution, and exhibited the system distribution in the form of a K-distribution. We underlined that this interesting distribution is highly relevant in some physical applications, and have emphasized some of its properties.

Fig. 1 .

 1 Fig. 1. Examples of K-distributions (12) with a = 0.1 to a = 10, and with b = 3a/2 chosen so that the expectation of all distributions is equal to 1.

  where Γ(a+1, 2 √ bx) is the upper incomplete Gamma function. It behaves as (2 √ bx) a e -2 √ bxfor large values of x and lim x→+∞ e λx F (x) = +∞ for all positive λ, which means that the distribution is heavytailed.

Fig. 2 .

 2 Fig.2. A typical realization of a K-distribution with a = 0.5.