Georges Da Costa
email: georges.da-costa@irit.fr

Ariel Oleksiak

Wojciech Piatek
email: piatek]@man.poznan.pl

Jaume Salom

Laura Sisó

Minimization of Costs and Energy Consumption in a Data Center by a Workload-based Capacity Management

Keywords: data centers, energy efficiency, simulations, heat-aware, metrics, OPEX, CAPEX

In this paper we present an approach to improve power and cooling capacity management in a data center by taking into account knowledge about applications and workloads. We apply power capping techniques and proper cooling infrastructure configuration to achieve savings in energy and costs. To estimate values of a total energy consumption and costs we simulate both IT software/hardware and cooling infrastructure at once using the CoolEmAll SVD Toolkit. We also investigated the use of power capping to adjust data center operation to variable power supply and pricing. By better adjusting cooling infrastructure to specific types of workloads, we were able to find a configuration which resulted in energy, OPEX and CAPEX savings in the range of 4-25%.

Introduction

The problem of capacity management in data centers is a well known issue, which data center planners and operators must deal with. The problem can be defined as finding such a data center configuration that its space, power and cooling capacity is maximized. In other words, the goal is to put maximal number of servers into a data center subject to its size, electrical infrastructure power limits, and heat dissipation constraints. Usually, this process is based on server power usage nameplates and by getting theoretical peak values from specifications. Unfortunately, these values are often the Power Supply Unit (PSU) maximum capacity so they substantially overestimate actual power loads. Therefore, vendors sometimes deliver calculators that help to obtain estimations closer to real values. Still, most of these methods neither take into consideration characteristics of specific applications nor dynamic properties of workloads that are executed in data centers. Some attempts to! apply more advanced power capping to improve efficiency of the whole data center can be found in literature. An alternative method to power capping based on managing distributed UPS energy is presented in [START_REF] Kontorinis | Managing distributed UPS energy for effective power capping in data centers[END_REF]. Interesting approach to combine IT workloads, power, cooling and renewable energy was studied in [START_REF] Liu | Renewable and cooling aware workload management for sustainable data centers[END_REF] but without use of power capping techniques. In [START_REF] Hankendi | vCap: Adaptive power capping for virtualized servers[END_REF] authors propose adaptive power capping for virtualized servers, however they investigate neither the cooling system nor variable power supply. Dynamic power capping to enable data center participation in power markets was proposed in [START_REF] Hao Chen | Dynamic server power capping for enabling data center participation in power markets[END_REF] but without detailed cooling consideration, either. To address these issues, we propose modeling and analysis of data center workloads and hardware to identify real power limits that should be met. Based on these limits we present methods to save energy and optimize cooling capacity of a data center including adaptation of limits to power supply and pricing. To meet this objective we have used the SVD Toolkit developed within the CoolEmAll project [START_REF] Dem Berge | Modeling and Simulation of Data Center Energy-Efficiency in CoolEmAll[END_REF]. The toolkit enables data center designers and operators to reduce its energy impact by combining the optimization of IT, cooling and workload management. For this purpose, CoolEmAll project investigated in a holistic approach how cooling, heat transfer, IT infrastructure, and applicationworkloads influence overall cooling-and energy-efficiency of data centers, taking aspects into account that traditionally have been considered separately. SVD Toolkit was used to conduct experiments described in this paper. In particular, most simulations were done using one of the main tools of the SVD Toolkit -the Data Center Workload and Resource Management Simulator (DCworms) [START_REF] Kurowski | DCworms -a tool for simulation of energy efficiency in distributed computing infrastructures[END_REF]. Using the CoolEmAll SVD Toolkit we demonstrate how to improve capacity management by taking into account knowledge about applications and workloads as well as by using power capping techniques and proper cooling infrastructure configuration. To obtain total energy consumption, we simulate both IT software/hardware and cooling infrastructure in parallel. In this way, by better adjusting cooling infrastructure to specific types of workloads, we were able to find a configuration which result in energy savings and even in improvement of CAPEX (Capital Expenditures) without significant workload performance deterioration. Decrease in CAPEX was achieved by the selection of smaller chiller which fits the foreseen workloads better. Energy savings were achieved by increase of server inlet temperature. This was possible by limiting power used by particular racks and by compliance to the latest ASHRAE recommendations. Finally, we applied power capping to adjust data center operation to ! variable power supply and achieved additional OPEX (Operating Expenditures) savings. The structure of this paper is as follows. In Section 2 we present a model of a data center including models of IT hardware, cooling, workloads and applications. This section also contains definitions of metrics used for the assessment of data center configurations studied in this paper. We analyze workloads along with their impact of on energy-efficiency in Section 3. Based on this analysis we define power limits which allow reducing energy consumption and costs of a data center operation. Section 4 contains results of the data center optimization using power capping methods and decisions about cooling infrastructure deployment and configuration. Section 5 concludes the paper.

Data center model

Modeling Workloads

In terms of workload management, workload items are defined as jobs that are submitted by users [START_REF] Feitelson | Workload Modeling for Computer Systems Performance Evaluation[END_REF]. Thus, modeling of workloads consists in providing information about structure, resource requirements, relationships and time intervals of jobs arriving to the management and scheduling system. Primary properties of a workload include:

number of jobs to be scheduled jobs arrival rate, expressed as a time interval between successive jobs reference to an application profile describing behaviour of particular job on the hardware (resource requirements and execution times)

The last one is described in the next section in more detail. Having these dependencies established, it is possible to express the impact of particular workload on the hardware layer. For now, one of the main and commonly used format that provides unitary description of workloads models and logs obtained from real systems is Standard Workload Format (SWF) [START_REF]ParallelWorkload Archive[END_REF].

As mentioned, workload profiles may be obtained by monitoring real systems or generated synthetically. The main aim of synthetic workloads is to reflect the behavior of real observed workloads and to characterize them at the desired level of detail. Moreover, they are also commonly adopted to evaluate the system performance for the modified or completely theoretical workload models. Usage of synthetic workloads and their comparison to the real ones have been the subject of research for many years [START_REF] Lo | A comparative study of real workload traces and synthetic workload models for parallel job scheduling[END_REF].

Modeling applications

Concerning application-led management a maximum feedback is needed from the applications from different point of view. The focus is on power-, energy-and thermal-impact of decisions on the system. Still it is impossible to put a wattmeter on an application. In order to obtain the same kind of information, we monitored applications to evaluate their resource consumption at each second. At each of these points, using system values and hardware performance counters, processor, memory and I/O resources are monitored. Using these information and models we produce for each of these timestamps an evaluation of the power consumption [START_REF] Da | Modeling the Energy Consumption of Distributed Applications[END_REF]. Each of those values are monitored, computed and stored in real-time in a database for future use.

In the system, an application is then described as the resources it uses on a particular hardware. Each application can be run on different hardware or configuration (frequency for example) and those data are associated with the same application. In case the data for a particular application on a particular hardware is not available, a translation tool is used to evaluate the behavior of the application using its behavior on a different hardware. First, it models the resource bottleneck of an application using the monitored resource consumption on a particular hardware. Using the target hardware specification, it evaluates the resource bottleneck and thus overall resource consumption on that hardware.

Using the monitored data, we create a description of applications based on their phases following the same methodology as in [START_REF] Landry | DNA-Inspired Scheme for Building the Energy Profile of HPC Systems[END_REF]. A phase is defined as a duration when resources consumption are stable. As an example, Fig. 1 show the profile of a Fast Fourier Transform algorithm with its phases. Using the XML files describing exact application behavior and resource consumption, SVD toolkit can evaluate precisely the impact of its decisions.

Fig. 1: Profile of the benchmark test3d: 3D real-to-complex FFT routine

Modeling servers

In the scope of CoolEmAll, data center server room is composed by a number of racks. Each rack consists of a set of node groups, which are then responsible for hosting a collection of nodes. Node groups are defined by a means of chassis that models the placement of nodes within the node group as well as mounted fans. The main component of the node is a processor with assigned number of cores and computing capability (expressed by a clock speed). Moreover, each processor comes with its power and computing profile, described by the means of C-States and P-States defining operating states with corresponding power usage values for different utilization levels. Node definition is supplemented by a description of memory and network. Rack represents a standardized enclosure for carrying server and power supply modules. Power profiles of IT infrastructure are the basis for calculating the power consumption of particular resources.

The following equations show how the power usage for different resource levels is estimated. P cpu (P x , load) = P cpu (P x , 0) + load * (P cpu (P x , 100) -P cpu (P x , 0))/100 (1)

where P cpu (P x , load) is a power consumed by a processor operating in a given P-State P x and utilized in a level denoted by load. P cpu (P x , 0) and P cpu (P x , 100) expresses an idle and fully loaded processor working in a given P-State, respectively (these constant values are part of the processor power profile providing power consumptions levels for all available frequencies).

P node = n i=1 P cpui + P mem + P net (2)
where P node is a power consumed by a node, n is the number of processors assigned to a node, P mem is a power drawn by a memory, while P net by a network.

P node group = m i=1 P nodei + k j=1 P f anj (3)
where P node group is a power consumed by a node group, m is the number of nodes placed in a node group, k is the number of fans mounted within it and P f anj is a power used by particular fan j.

P rack = ( l i=1 P node groupi )/η psu (4)
where P rack is a power consumed by a rack, l defines the number of carried node groups and η psu is efficiency of a power supply unit.

Finally, each component is accompanied with its carbon emissions and electricity costs. Apart from IT equipment, data center server room is composed by a cooling devices, which are the subject of next subsection.

Cooling models

The SVD CoolEmAll toolkit integrates models to calculate the power associated to cooling equipment and other electric facilities required in data center to fulfill its mission related with IT services. The cooling model provided consists of a simple data center where central fan and air-water coil cools the IT equipment and other related loads (PDU, UPS and lighting). A chiller placed outside provides cooling water to the coil and dissipates the exhausted heat from the room to the atmosphere by a dry-cooler (Figure 2 shows details). The power model adds the consumption of IT, fans, chiller, PDU and lighting. Other electric components of a data center as back-up generator or transformer are excluded from the present model.

The following model description is based on a single time-stamp where Q is referred to heat dissipated and P to power consumption. The time variability is indicated by (t). This model has been constructed based on basic thermodynamic equations of conservation of mass and energy. The total power consumption of a data center (P DC ) will be calculated with Eq. 5, where P load DC is the power used by IT components, P chiller is the consumption of the chiller, P f ans DC is Fig. 2: Model of cooling and power facilities of a data center the consumption of fans in data center and P others is the consumption associated to PDU and lighting:

P DC (t) = P load DC (t) + P chiller (t) + P f ans DC (t) + P others (t) (5)
The total thermal load (Q DC ) is the sum of the heat associated to IT load (Q load DC ), the heat from other loads, as PDU and lighting (Q others DC ) and the heat from fans distributing air inside a data center room (Q f an DC ).

Q DC (t) = Q load DC (t) + Q f an DC (t) + Q others DC (t) (6)
The cooling demand that should be covered (Q cooling ) is the thermal load in data center including the inefficiencies in the air-water coil represented by η cc according Eq. 7. That corresponds to the heat exchanger efficiency of a common CRAH, where heat of the room is transferred to the water flow (Q cooling ).

Q cooling (t) = Q DC (t) η cc (7)
The chiller has been modelled with generic profiles based on condenser temperature (T co ), evaporator temperature (T ev ) and partial load ratio (P LR). Thereby, the model presented here should provide a general method to determine the power consumption of the chiller without knowing the specific characteristics of the chiller provided by a certain manufacturer. As a result, it has been used parametric curves implemented by the Building Certification Code in Spain [1] named COOL(T ev , T co ) and CoolP R(T ev , T co , P LR, EER rated ) following certain relations depicted in Eq. 8 and Eq. 9.

Q cooling nom = Q cooling rated • COOL(T ev , T co ) (8) CoolP R(t) = CoolP R(T ev , T co , P LR, EER rated ) = 1 EER(t) (9)
Partial Load Ratio is the relation between the cooling demand in a certain conditions and the cooling load in nominal conditions (Q cooling nom ) corresponding to the operation of the chiller at the chilled water temperature (T ev ) and condenser water temperature (T co ) set-up (Eq. 10). At the same time, Q cooling nom has relation with the cooling capacity rated (Q cooling rated ) which corresponds to load of the chiller in Standard Conditions (full load; temperature of chilled water leaving the chiller at 7 o C and temperature of condenser water entering the chiller at 30 o C) as stated in Eq. 8.

P LR(t) = Q cooling (t) Q cooling nom (10)
The the relation between the cooling load and the power consumed in the chiller (P chiller ) is linked by the Energy Efficiency Ratio (EER), that quantifies the cooling provided by the chiller by each unit of power consumed, according Eq. 11. EER rated corresponds to the value of the parameter measured at Standard Conditions defined above.

P chiller (t) = Q cooling (t) EER(t) (11)

Assessment of data center efficiency, performance, and costs

Metrics CoolEmAll SVD Toolkit provides a set of metrics divided in the level of granularity of the analysis (node, node-group, rack and data center). The whole group of metrics assesses the resource usage, capacity, energy, heat-aware, green and financial concepts. The total selection of metrics of CoolEmAll are described in public report of the project [START_REF] Sisó | D5.1 White paper on Energy-and Heataware metrics for computing modules -CoolEmAll Deliverable[END_REF] as well as in some articles [16] [17]. Total Energy Consumed: this corresponds to the total energy consumed by the data center in a certain period of time.

Power Usage Effectiveness (PUE): defined by The Green Grid [START_REF] Avelar | The Green Grid[END_REF] this metric consist of dividing power used by the data center between power used by the IT equipment. The accuracy level of the metrics is related with the point of measurement of IT power, that can be the UPS (Uninterruptible Power Supply Unit), the PDU (Power Distribution Unit) or the IT itself, after PSU (Power Supply Unit). When the measurement is done after the PSU the metric is defined as PUE Level 3.

When the measurement is referred to IT properly, excluding PSU and fans, the metric is named PUE Level 4, according CoolEmAll project proposal [START_REF] Sisó | Energy and Heat-Aware Metrics for Data Centers: Metrics Analysis in the Framework of CoolEmAll Project[END_REF].

Carbon emissions: this metric is calculated multiplying the total power consumed by carbon emissions factor (CEF). CEF depends on the country power generation mix and power system efficiency. For the approach of this study, 0.34kg/kWh has been used as average value for the European Union according to [START_REF] Kemma | Methodology Study Eco-design of Energyusing Products MEEUP[END_REF] OPEX: it is calculated multiplying the total power consumed by the price of electricity. The price of electricity has been considered as 0.0942 e/kWh for EU-28 as average of 2013 according to [START_REF]European Commission[END_REF].

CAPEX: it is the amount of money used to acquire equipment or to improve the useful life of existing facilities.

Analysis of workloads

As mentioned in Section 2.1, workloads are characterized by the number of jobs, their arrival rate, resource requirements and execution time of particular applications. The following section contain describes the results of workload simulations performed by the means of Data Center Workload and Resource Management Simulator, which is part of SVD Toolkit.

Simulation of diverse workloads using DCworms

Resource characteristics In our experiments we used a configuration of the real server room. Each server was equipped with a processor belonging to Intel Xeon processors family. The following table (

Identifying power caps

Based on the simulation results obtained for execution of both workloads using Load Balancing policy, we observed two visible increases on utilization criteria, reaching almost 75% and 95% in the highest peak for Workload 1 and Workload 2 respectively. High utilization values have direct impact on the power consumption and thus might result in sudden power drawn peaks. Identification of such levels is crucial in terms of avoiding hot spots and decreasing data center costs. Taking into account power consumption ranges for the modeled server room, power consumption distribution obtained during the experiments and the utilization curves we decided to use the following approach to specify the values of power caps. As there occured temporary, but significant load rises and we were not considering the possibility of switching nodes on/off, we wanted to ensure constant computational capabilities for all the servers within particular racks. To this end the power cap level is determined by the total power consumption of the rack, with all the processors fully loaded and working in the highest P-State (with lowest frequency). The following formula can be used to calculate this value (P C) for the given rack j.

P C j = n i=1 P CP Ui (P hi , 100%), ( 12
)
where n is the number of processors in a rack, P CP U is the power consumed by the processor working under given utilization level and in the given P-State, P h refers to the highest P-State (power consumption is lower at higher P-State).

On the other hand, in order not to observe the performance losses (due to frequency downgrading) another threshold is necessary. It aims at setting the power consumption level P U below which the current processor performance state will increase. It is defined by the following equation:

P U j = P C j • n i=1 P CP Ui (P hi , 100%) n i=1 P CP Ui (P hi-1 , 100%) , ( 13
)
where n is the number of processors in a rack, P CP U is the power consumed by the processor working under given utilization level and with the given P-State, P h and P h-1 refer to two highest P-States. As power consumption is lower at higher P-States, thus, P U j is lower than P C j .

P U j allows increasing the current processors performance states at least by one without exceeding the power cap limit (P C j ). Below Adjusting power limits to workloads Having information about historical or predicted workloads it is possible to adjust power caps. For instance, there may exist specific patterns of incoming tasks related to peak hours, time of a day, etc. This knowledge can be applied to identification of optimal power caps.

There are two main requirements that should be taken into consideration while setting the values of power caps. First of all, the use of power capping shouldn't cause significant increase of IT energy consumption for a given workload. Second, the mean completion time of tasks should not go below certain required threshold.

The first requirement can be defined as follows. Let assume that energy decrease caused by power capping in rack j is denoted as E excess j and given in Eq. 14. This amount of energy can be illustrated by the field above the power cap line in Figure 3. On the other hand, let denote by E reserve j the amount of additional energy that can be used in a rack without exceeding the set power cap. This can be seen as a free space below the power cap in Figure 3 and defined by Eq. 15.

E excess j = t2 t1 max(0, P IT j (t) -P C j )dt ( 14
)
E reserve j = t2 t1 max(0, P C j -P IT j (t))dt ( 15
)
Fig. 3: Power distribution for Workload 1 and Workload 2 for two racks Then, the condition E excess j < E reserve j must be met. Otherwise, tasks whose execution times are increased by decreased performance states of CPUs could cause additional delays of additional tasks. Of course, this method is approximation as the actual results depend on sizes of tasks, their distribution, and exact relation between CPU performance states and execution time. However, as results in Section 4.3 show this approach is helpful to avoid increase of energy consumption caused by power capping.

The second requirement is meant to limit the power but without visible performance lost. The mean completion time increase caused by power capping can be estimated as a product of the CPU frequency change (we assume proportional relation to execution time) and a percentage of time for which power capping was used. We empirically set a 5% as a threshold for mean completion time increase to limit overall delays of the workload completion time and this condition was met (see results in Table 6). This parameter was used to limit CPU frequency decrease according to a model presented in Section 2.3 and can be based on specific Service Level Agreements with end users. [START_REF] Hao Chen | Dynamic server power capping for enabling data center participation in power markets[END_REF] Optimizing capacity using power capping methods

Power capping methods

Generally, power capping solutions can be divided into: software-based (coarsegrained and slower) and hardware-based solutions (fine-grained and faster).

Software-based solutions can be introduced independently from the vendor and regardless of whether hardware power capping is available. It can be applied on higher levels, e.g. managing tasks in a queue and balancing the load (with respect to power) among racks. The drawback of the software-based approach is longer time of reaction and more coarse-grained granularity.

Hardware-based power capping addresses this issue by the means of two main technologies available at processor level that enable the use of power capping. The first one is related to processor P-States and consists in lowering the processor core frequency and voltage. That provides a good power reduction for a relatively small loss in performance. However, using P-States can lower power consumption only to a certain point. Reducing consumption below that point requires the use of second technology, namely clock throttling. In this case, depending on the processor model, the system BIOS can either reprogram the processor to run at a lower frequency or modulate the processor between running periods and stopped periods.

In this paper we focus on the hardware-based approach benefiting from the processors P-States, as in the real data centers it ensures more reliable and faster effects. Moreover, it is often supported by hardware vendors and can be easily applied on the resource management level without affecting existing queueing system configuration (comparing to software-based approach). Its pseudo code for a rack is depicted by Algorithm 1.

Simulation experiments

To study the impact of power capping approach we performed another two simulations each time applying the power caps levels introduced in Section 3.2. Moreover, in the first simulation run we increased the inlet temperature (temperature of air entering the room) to 27 o C, while in the latter one we additionally modified, according to Table 5, the cooling capacity rated factor to 180 [kW ]. Below we introduce the nomenclature used to compare the simulation results.

Simulation results

This section shows the simulation results for three types of experiments performed. Due to the paper constraints, only the results for Workload 1 are presented (Table 6). Total savings achieved The following section shows the savings achieved for the simulations carried out. The reference case, with Load Balancing policy, is named "A". Optimized cases are named "B" and "C" respectively. First, when the strategy of power capping is applied, main savings can be observed in the chiller consumption (reaching 37,50% and 17.09% respectively) due to the efficiency of the chiller (EER) improves with higher inlet temperatures. This leads to savings in terms of total energy consumed that are equal to 4.19% in case "B'" and 1.20% in case "C".

In the simulation carried out, the strategy consisted of cutting the maximum power of racks keeping same cooling facilities (in case "B") or changing the chiller capabilities (in case "C"). The result obtained in these cases is a reduction in OPEX associated to power saved mainly in chiller and in CAPEX due to reduction of IT infrastructure. The metrics calculated from the results of those simulations are shown in Table 7. Proposed approaches provide small benefits on PUE -the savings are obtained due to the lower power consumption of the chiller. With a power capping of 10% the savings obtained in total energy consumed, carbon emissions and electricity costs (OPEX) are 4.30% for case (B) and 1.29% for case (C). The corresponding values obtained in savings extrapolated to a whole year considering a 24x7 operation time are 60MWh/year, 20 tones CO 2 /year, 5666 Euros/year and 21MWh/year, 7 tones CO 2 /year, 1982 Euros/year, respectively. Also, the CAPEX costs associated to less equipment required are calculated based on the following approach. Total building cost of traditional data center is estimated as 15 million-US$ per MW of IT load according market survey developed by 451 Research company, referred as [19]. Converting this value to Euros with average annual ratios determined by the European Central Bank [21] referred to 2012, the corresponding value is 10784 Euros/kW(IT). On the other hand, the following distribution of cost between subsystems is considered according the study done by Schneider Electric [START_REF] Rassmussen | Determining Total Cost of ownership for data center and network room infrastructure[END_REF].

Metrics

The 10% capping on maximum power of racks will affect directly the cost of those IT equipment but also on the sub-system of power equipment. Table 8 shows the distribution of costs of the three cases simulated. The costs of case (B) and (C) have been calculated estimating a reduction of 10% in racks and power equipment. Finally, with this assumption, the savings obtained in CAPEX over the total cost of the data center is a 4% or 62 thousands of Euros and 7% or 109 thousands of Euros, respectively.

Application to demand-response management

Nowadays, power grids face significant transformations. More open energy market, increased contribution of renewable energy sources, and rising energy prices stimulate changes of power grids to cope with new challenges such as adaptation to changing demand and supply, i.e. demand-response management. The approach to apply demand-response management to data centers was also already studied, e.g. proposed in [START_REF]The All4Green project website[END_REF]. We show that our approach to analysis of workloads and power capping mechanism can be applied to reduce costs in data centers.

Let's assume that for a period assumed in previous Sections (3h 20min) there is a regular price for energy: 0.0942/kWh. Now, let's also assume that period of the same size is a peak period in which energy provider is struggling with a demand that exceeds provider's supply. The provider to cope with this demand proposes the following contract to its customers: a regular price for this period will stay on the same level provided that a customer guarantees that it will not exceed 200kW of power at anytime. Otherwise, the cost of 1 kWh will rise up to 0.15/kWh. To reduce costs in this case we applied power capping to the peak period. The comparison of approaches without and with power capping are presented in Table 9. In the first case power capping was not used in any period. In the second case power capping was applied to the second (peak) period. As it can be easily seen, the total cost savings reached almost 25%. Extrapolating these numbers to the whole year would give around 45000e of savings.

Approaches

Conclusions

In this paper we demonstrated the use of the CoolEmAll SVD Toolkit to improve power and cooling capacity management in a data center by taking into account knowledge about applications and workloads. We applied power capping techniques and proper cooling infrastructure configuration to achieve savings in energy and costs. To obtain estimated values of a total energy consumption we simulated both IT software/hardware and cooling infrastructure using our tools. In this way, by better adjusting cooling infrastructure to specific types of workloads, we were able to find a configuration which resulted in energy savings by around 5% and corresponding OPEX decrease. We have also found improvements of CAPEX without significant workload performance deterioration. Decrease in CAPEX was achieved by the selection of smaller chiller which is sufficient for the foreseen types of workloads. Savings in CAPEX reached 7% for the case in which a smaller chiller was used according to the work! load analysis results and power capping strategies. Replacing only electrical equipment brought 4% of savings in CAPEX. Energy savings were achieved by increase of the server inlet temperature. This was possible by limiting power used by particular racks and by compliance to the latest ASHRAE recommendations. Finally, we applied power capping to adjust data center operation to variable power supply and pricing. We achieved additional OPEX savings in order of 25% (45000eper year in the studied case).

Future work will include further improvements and tuning of cooling models. It will also include closer integration of CFD simulations into this analysis in order to identify hot spots and other consequences of modifications in a data center configuration. This approach will be used for various types of data centers. Finally, we plan to study more dynamic power capping strategies by adjusting power caps to the situation in a data center such as level and priority of load, energy supply and prices.

-

  Experiment A: Load Balancing strategy, T R in = 18 o C, reference case. -Experiment B: Load Balancing with Power Capping approach, T R in = 27 o C -Experiment C: Load Balancing with Power Capping approach, T R in = 27 o C, Q cooling rated = 180[kW ]

Fig. 4 :

 4 Fig. 4: Example power distribution on two racks for Workload 1 before (left) and after (right) applying a power capping technique

Table 1 )

 1 summarizes overall characteristics of particular racks.

	Rack name Number	Number of	Processor type Min. power us-	Max. power usage
		of nodes	processors		age (idle) [W]	(100% load) [W]
	Rack 1	84	2	Xeon E5-2603 10292	27672
	Rack 2	84	2	Xeon E5-2630 12030	30568
	Rack 3	84	1	Xeon L5310	4499	11258
	Rack 4	84	1	Xeon L5310	4499	11258
	Rack 5	84	2	Xeon E5-2603 10292	27672
	Rack 6	84	2	Xeon E5-2603 10292	27672
	Rack 7	84	2	Xeon E5-2630 12030	30568
	Rack 8	56	2	Xeon E5-2630 8020	20379
	sum	644	1120	-	71955	187046

Table 1 :

 1 Power characteristics of racks in the server room Additionally, server room was equipped with the cooling facilities presented in Table2. Temperature difference between Tev and TR in ∆T hex 10 o C

	Parameter	Symbol in the	Value
		equations
	Cooling capacity rated	Q cooling rated 240000 [W]
	Energy efficiency ration rated	EER rated	3
	Efficiency of cooling coil	ηcc	0.95
	Data center fans efficiency	η f	0.6

Table 2 :

 2 Cooling facilities characteristicsFinally, the following input parameters were applied to the simulation environment (Table3).

	Parameter	Symbol in the	Value
		equations
	Relation between P loadDC and P others α	0.2
	Inlet temperature	TR in	18 o C
	Outlet temperature	TR out	33 o C
	Pressure drop	∆p	65 J/m 3

Table 3 :

 3 Input parametersWorkloads and application profiles In our experiment we evaluated two workloads with different utilization levels what was achieved by the modification of arrival rate (all tasks arrive according to the Poisson distribution) and the number of submitted tasks. The former workload consists of 1280 tasks, while the latter consists of 1760 tasks.A distribution of applications constituting both workloads is the same in both cases and looks as follows: 20% -App1, 50% -App2, 30% -App3. Their general overview is shown in Table4. The understanding of the cells content is as follows: number of requested processors, execution time, load level (in [%]).

	Processor Type	App1	App2	App3
	Xeon E5-2630	1, 380, 84	4, 3200, 62.6	6, 3200, 94
	Xeon E5-2603	1, 400, 86	4, 3600, 92	-
	Xeon L5310	1, 1200, 92	-	-

Table 4 :

 4 Application characteristics

  table introduces boundary values according to the aforementioned approach;

	Rack name Rack 1 Rack 2 Rack 3 Rack 4 Rack 5 Rack 6 Rack 7 Rack 8 sum
	PC level [W] 20333 21878 8940	8940	20333 20333 21878 14585 137220
	PU level [W] 18854 20223 8449	8449	18854 18854 20223 13482 127388

Table 5 :

 5 Power caps values for the racks in the server room

  Algorithm 1 Pseudo code of power capping algorithm State=higher frequency for each processor i in a rack with Px i = Py do Px i = Px i +1 if P <= P C then break end if end for until P <= P C else if P < P U then repeat Py ← highest P-State of all processors ⊲ higher P-State=lower frequency for each processor i in a rack with Px i = Py do Px i = Px i -1 if P >= P C then break end if end for until P >= P C or Py= lowest available P-State of all processors end if

	Require: P	⊲ description of current power consumption of a rack
	Require: P C		⊲ power cap level for a rack
	Require: P U		⊲ power threshold for a rack
	Ensure: Px i		⊲ final P-state of a processor i
	if P > P C then
	repeat
	Py ← lowest P-State of all processors	⊲ lower P-

Table 6 :

 6 Simulation results for Workload 1 Figure4depicts the power distribution before and after applying a power capping technique.

Table 7 :

 7 Savings on particular metrics

Table 8 :

 8 Cost of data center placed in a room for three cases (thousand Euro)

	Costs by sub-system	A	B	C
	project management	156	156	156
	power equipment	562	506	506
	cooling equipment	187	187	141
	engineering & installation 562	562	562
	racks	62	56	56
	system monitoring	31	31	31
	TOTAL	1562	1500	1453

Table 9 :

 9 Comparison of approaches with and without power capping to deal with high demand periods

		Total	energy	Average	energy	Mean completion
		cost [e]		price [e]		time [s]
	no power capping 128.24		0.12		6919
	mix	96.8		0.0942		7090

Acknowledgements. The results presented in this paper are partially funded by the European Commission under contract 288701 through the project CoolE-mAll and by a grant from Polish National Science Center under award number 2013/08/A/ST6/00296.


Alin Achim

Adrian Basarab

George Tzagkarakis

Panagiotis Tsakalides

Denis Kouamé

RECONSTRUCTION OF COMPRESSIVELY SAMPLED ULTRASOUND IMAGES USING DUAL PRIOR INFORMATION

Keywords: Compressive sensing, ℓ p norm, alphastable distributions, medical ultrasound

This paper introduces a new technique for compressive sampling reconstruction of biomedical ultrasound images that exploits two types of prior information. On the one hand, our proposed approach is based on the observation that ultrasound RF echoes are best characterised statistically using alpha-stable distributions. On the other hand, through knowledge of the acquisition process, the support of the RF echoes in the Fourier domain can be easily inferred. Together, these two facts inform an iteratively reweighted least squares (IRLS) algorithm, which is shown to outperform previously proposed reconstruction techniques, both visually and in terms of two objective evaluation measures.

INTRODUCTION

Recent developments in medical ultrasound (US) imaging have led to commercial systems with the capability of acquiring Real-Time 3D (RT3D or 4D) image data sets. However, current scanners can only produce a few volume images per second, which is fast enough to see a fetus smile but not fast enough to see heart valves moving. Compressive sensing (CS) could prove to be a powerful solution to enhance US images frame rate by decreasing the amount of acquired data. In this context, in the last four years, a few research groups worked on the feasibility of compressive sampling in US imaging and several attempts of applying the CS theory may be found in the recent literature (for an overview see e.g. [START_REF] Liebgott | Compressive sensing in medical ultrasound[END_REF]).

In particular, in [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF], we introduced a novel framework for CS of biomedical ultrasonic signals based on modeling data with symmetric alpha-stable (SαS ) distributions. Then, we proposed an ℓ p norm-based minimization approach that employed the iteratively reweighted least squares (IRLS) algorithm, but in which the parameter p was judiciously chosen by relating it to the characteristic exponent of the underlying This work was partially funded by the CS-ORION (PIAP-GA-2009-251605) grant within the 7th Framework Program of the European Community alpha-stable distributed data. The results showed a significant increase of the reconstruction quality when compared with previous ℓ 1 minimization algorithms. On the other hand, the effect of the random sampling pattern on the reconstruction quality, when working in the frequency domain (k-space) was studied in [START_REF] Quinsac | Compressed sensing of ultrasound images: Sampling of spatial and frequency domains[END_REF]. This was further exploited in [START_REF] Basarab | Medical ultrasound image reconstruction using compressive sampling and l p norm minimization[END_REF] for the design of a US reconstruction technique similar to [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF] but operating in the Fourier domain.

In this contribution, we further extend our techniques described in [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF][START_REF] Basarab | Medical ultrasound image reconstruction using compressive sampling and l p norm minimization[END_REF] by supplementing the prior information available to an ℓ p norm minimisation algorithm with the support of the RF echoes in the frequency domain. In ultrasound applications the latter can be easily inferred through knowledge of the acquisition frequency and transducer bandwidth.

In the following section, we provide a brief, necessary overview of the compressive sensing theory and of the heavytailed model that we employ for ultrasound data. Section 3 describes our proposed method that exploits dual prior information, while Section 4 illustrates its reconstruction performance. Finally, Section 5 summarises the main results.

BACKGROUND

Compressive sensing

It has been demonstrated [START_REF] Candes | Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information[END_REF][START_REF] Donoho | Compressed sensing[END_REF] that if a signal is K-sparse in one basis (meaning that it can be represented by K elements of that basis), then it can be recovered from M = c.K << N fixed (non-adaptive) linear projections onto a second basis, called the measurement basis, which is incoherent with the sparsity basis, and where c > 1 is a small overmeasuring constant. The measurement model is

y = Φx (1)
where x is the N × 1 discrete-time signal, y is the M × 1 vector containing the compressive measurements, and Φ is the M × N measurement matrix. Using the M measurements in the first basis and given the K-sparsity property in the other basis, the original signal can be recovered by taking a number of different approaches. In the context of this work, our interest lies in non-convex optimisation (re-weighted ℓ p minimisation [START_REF] Chartrand | Iteratively reweighted algorithms for compressive sensing[END_REF][START_REF] Miosso | Compressive sensing reconstruction with prior information by iteratively reweighted least squares[END_REF]) strategies.

α-stable distributions as models for RF echoes

The ultrasound image formation theory has been long time dominated by the assumption of Gaussianity for the return RF echoes. However, the authors in [START_REF] Kutay | On modeling biomedical ultrasound RF echoes using a power-law shot noise model[END_REF] have shown that ultrasound RF echoes can be accurately modelled using a powerlaw shot noise model, which in [START_REF] Petropulu | Power-law shot noise and its relationship to long-memory alphastable processes[END_REF] has been in turn shown to be related to alpha-stable distributions.

By definition, a random variable is called symmetric αstable (SαS) if its characteristic function is of the form:

ϕ(ω) = exp( jδ ω -γ|ω| α ), (2)
where α is the characteristic exponent, taking values 0

< α ≤ 2, δ (-∞ < δ < ∞)
is the location parameter, and γ (γ > 0) is the dispersion of the distribution. For values of α in the interval (1, 2], the location parameter δ corresponds to the mean of the SαS distribution, while for 0 < α ≤ 1, δ corresponds to its median. The dispersion parameter γ determines the spread of the distribution around its location parameter δ , similar to the variance of the Gaussian distribution.

The alpha-stable tail power law provided one of the earliest approaches in estimating the stability index of real measurements [START_REF] Samorodnitsky | Stable Non-Gaussian Random Processes: Stochastic Models with Infinite Variance[END_REF]. The empirical distribution of the data, plotted on a log-log scale, should approach a straight line with slope α if the data is stable. Maximum likelihood methods developed by various authors are asymptotically efficient and have become amenable to fast implementations [START_REF] Nolan | Maximum likelihood estimation of stable parameters[END_REF]. In Fig. 1 we show an example of an ultrasound RF echo modelled using SαS density functions both in time and in frequency domain.

IRLS WITH DUAL PRIOR INFORMATION

In reconstructing an ultrasound RF echo, we would ideally like to be able to reconstruct a vector x, with the smallest number of non-zero components, that is,with the smallest ℓ 0 norm. Since this is however an NP-hard problem, several suboptimal strategies have been proposed with the majority solving a constrained optimization problem by employing ℓ 2 or ℓ 1 norms. CS reconstruction methods were also developed (e.g. [START_REF] Chartrand | Iteratively reweighted algorithms for compressive sensing[END_REF][START_REF] Miosso | Compressive sensing reconstruction with prior information by iteratively reweighted least squares[END_REF]) by employing ℓ p norms with p < 1, with the goal of approximating the ideal ℓ 0 case. Specifically, the problem consists in finding the vector x with the minimum ℓ p norm by minimising x = min x p subject to Φx = y

(3)

In [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF] we have devised a principled strategy for choosing the optimal p by relating the ℓ p norm minimisation to the actual SαS statistics of the RF signals. We achieved that by observing that for alpha-stable signals, which do not possess finite second-or higher-order moments, the minimum dispersion criterion [START_REF] Samorodnitsky | Stable Non-Gaussian Random Processes: Stochastic Models with Infinite Variance[END_REF] can be defined as an alternative to the classical minimum mean square error for Gaussian signals. This leads to a least ℓ p norm estimation problem, an approach that we have shown to enhance the reconstruction of heavy-tailed RF signals from their measurement projections [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF].

Here, our approach to RF signal reconstruction still relies on SαS-IRLS [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF] but is implemented in the frequency domain as in [START_REF] Basarab | Medical ultrasound image reconstruction using compressive sampling and l p norm minimization[END_REF] and modified (following [START_REF] Miosso | Compressive sensing reconstruction with prior information by iteratively reweighted least squares[END_REF]) to incorporate information on the support of RF signals.

Specifically, denote by X ∈ R N×J an US RF image formed by J RF signals of length N, x 1 , x 2 , . . . x J . The 1D Fourier transforms of all individual RF echoes x j can be written as

ξ j = F x j , j ∈ {1, 2, . . . J} (4)
where F ∈ C N×N is the 1D Fourier matrix. In the frequency domain, the measurement model becomes

m j = Φ j ξ j = Φ j F x j , j ∈ {1, 2, . . . , J} (5)
where Φ j are Gaussian matrices of size M × N (M ≪ N) and m j ∈ C M×1 . Now denote by Θ j the subset of points in {1, 2, . . . N} that defines the support of ξ j :

ξ j,k = 0, ∀k ∈ Θ j , j ∈ {1, 2, . . . , J} (6)
Following the arguments in [START_REF] Miosso | Compressive sensing reconstruction with prior information by iteratively reweighted least squares[END_REF], the information represented by ( 6) can be added to the IRLS algorithm for the minimisation of the l p norm by solving the following problem instead of (3) (or its frequency domain equivalent)

min ξ 1 2 N ∑ k=1 k / ∈Θ ξ j,k p subject to Φ j ξ j = m j , j ∈ {1, 2, . . . , J} (7)
Intuitively, (7) will offer a better solution than (3) because it will attempt to minimise the number of nonzero elements in ξ only outside the set Θ.

To solve [START_REF] Chartrand | Iteratively reweighted algorithms for compressive sensing[END_REF] we use the modified IRLS algorithm proposed in [START_REF] Miosso | Compressive sensing reconstruction with prior information by iteratively reweighted least squares[END_REF]. Following the arguments in [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF], the parameter p is set as α -0.01, where α is obtained by fitting an alpha-stable distribution to the data. Finally, the reconstructed RF lines are obtained by inverting the corresponding Fourier transforms:

x j = F -1 ξ j , j ∈ {1, 2, . . . , J} (8)

RECONSTRUCTION RESULTS

In this section we present reconstruction experiments conducted using real data corresponding to an in vivo healthy thyroid. The image was acquired with a Siemens Sonoline Elegra scanner using a 7.5 MHz linear probe and a sampling frequency of 50 MHz. Various sections of the original image were cropped and patches of size 256 × 512 were obtained. These patches were then sampled line by line using linear projections of random Gaussian bases at two levels. The two levels correspond to the number of samples taken from the original signal (the echo lines); these are 33% and 50% (i.e. M = 0.33N and M = 0.5N). Reconstruction of the samples was achieved by using the proposed ℓ p -norm minimization scheme (as described in Section 3) and for comparison, reconstructions using ℓ pnorm minimization with SαS -IRLS [START_REF] Achim | Compressive sensing for ultrasound RF echoes using α-stable distributions[END_REF] and SαS-IRLS in the Fourier domain (FD-SαS-IRLS) [START_REF] Basarab | Medical ultrasound image reconstruction using compressive sampling and l p norm minimization[END_REF] are shown (Table I). The values of α for each line and so that of p (which is derived from α), were estimated directly from the ultrasound RF signal while for our new approach (IRLS-DP) the support was inferred through knowledge of the frequency of acquisition and transducer bandpass as detailed above.

An analysis of the results was undertaken in terms of reconstruction quality, which was measured by means of the structural similarity index (SSIM) [START_REF] Wang | Image quality assessment: From error visibility to structural similarity[END_REF] and normalised root mean squared error (NRMSE) of the reconstructed echoes ensemble compared with the original ensemble. SSIM resembles more closely the human visual perception, and as such, it is often preferred than the commonly used MSE-based metrics. For a given image I and its reconstruction Î the SSIM is defined by: SSIM = (2µ

I µ Î + c 1 )(2σ I Î + c 2 ) (µ 2 I + µ 2 Î + c 1 )(σ 2 I + σ 2 Î + c 2 ) , (9)
where µ I , σ I are the mean and standard deviation of I (simi- larly for Î), σ I Î denotes the correlation coefficient of the two images, and c 1 , c 2 stabilize the division with a weak denominator. In particular, when SSIM equals 0 the two images are completely distinct, while when the two images are matched perfectly SSIM is equal to 1.

It can be seen in Table I that according to both metrics employed, the best results are obtained using the proposed reconstruction algorithm, which exploits two types of prior information. The results support the fact that reconstructing ultrasound RF echoes in the Fourier domain produces better results than directly in time domain. We attribute this to the sparser representations that can be achieved for RF lines in the Fourier domain. Taking into account prior information of the signal in the form of its support is also confirmed to optimise reconstruction. We should note however that, unlike the observation made in [START_REF] Miosso | Compressive sensing reconstruction with prior information by iteratively reweighted least squares[END_REF], reducing further the order p leads actually to worse reconstruction results. The optimal value of p is thus confirmed to be related to the underlying alphastable distribution of the data.

For a qualitative analysis, reconstruction results obtained with the three schemes discussed in this paper, with M = 0.33N measurements are also presented in Fig. 2. Visually, it can be seen that the IRLS-DP reconstruction introduces the least distortion, clearly producing the best result compared to the original and confirming the results indicated by the NMRSE and SSIM values obtained.

CONCLUSIONS

In this paper, we extended our previously proposed framework for compressive ultrasound imaging. We have shown through simulations that RF echoes can be best reconstructed by driving an ℓ p -norm minimisation problem with dual prior information: the value of the characteristic exponent of the RF line and its sparse support in the frequency domain. The latter plays certainly a central role in allowing a significant reduction in both the number of required measurements and computational cost. Nevertheless, our experiments strongly suggest that the optimal value of p in a ℓ p -norm minimisation procedure shouldn't be arbitrarily small but rather close to the characteristic exponent of the underlying alpha-stable distribution of the data. Our results are general but have been illustrated through only one example dataset due to space limitations.

Fig. 1 .

 1 Fig. 1. Example RF signal modelling with SαS distributions. (a) An RF signal in time domain. (b) The real part of its 1D Fourier transform.

Fig. 3 .

 3 Fig.3. Reconstruction errors for one RF line sampled at 33%. Top to bottom: original signal and reconstructed using SαS-IRLS, SαS-IRLS in the Fourier domain, and IRLS with dual prior respectively. Left column: RF lines; right column: the corresponding errors.

Fig. 2 .

 2 Fig. 2. Reconstruction results for a thyroid ultrasound image using 33% of the number of samples in the original. (a) B-mode ultrasound image. (b) SαS-IRLS reconstruction. (c) SαS-IRLS in the Fourier domain. (d) Fourier domain IRLS with dual prior.

Table 1 .

 1 Objective evaluation of three IRLS-based methods for reconstructing ultrasound images from RF frames with sampling rates of 33% and 50% relative to the original.

	M N %	Metric	Method SαS-IRLS FD-SαS-IRLS IRLS-DP
	33	NMSRE SSIM	0.697 0.208	0.540 0.586	0.249 0.908
	50	NMSRE SSIM	0.518 0.377	0.291 0.844	0.158 0.944


Jessica Bertheloot

François Barbier

Frédéric Boudon

Maria-Dolores Perez- Garcia

Thomas Péron

Christine Beveridge

Christophe Godin

Soulaiman Sakr

Sucrose interacts with auxin in the burst of axillary buds

Research focus


est destinée au dépôt

Modelling physiological regulation of bud burst by sucrose and auxin

Jessica Bertheloot, François Barbier, Frédéric Boudon, Maria-Dolores Perez-Garcia, Thomas Péron, Christine Beveridge, Christophe Godin,

Soulaiman Sakr

To cite this version:

Jessica Bertheloot, François Barbier, Frédéric Boudon, Maria-Dolores Perez-Garcia, Thomas Branching is difficult to control in production systems due to the lack of knowledge on the mechanisms regulating branching plasticity and response to environment. Our research aims at modelling the physiological mechanisms responsible for the spatio-temporal patterns of bud burst using rosebush as a plant model.

The burst of an axillary bud is inhibited by auxin, and several decades of research have led to the identification of the physiological pathway by which auxin regulates bud burst (Rameau et al., 2015). Recent studies have enlightened the importance of a novel player, sugars, and that sucrose modulates hormonal mechanisms suggesting a non-trivial interaction between sucrose and auxin (Mason et al., 2014;Barbier et al., 2015). The present study aims at modelling how sucrose interacts with auxin and other key hormones to regulate the burst of a bud. The "bud" model developed will serve as an elementary brick for a FSPM model simulating bud burst patterns on a plant.

Methods.

Current knowledge allows defining a model of auxin action on the bud: auxin level in the stem controls the level of two "second messenger" hormones, cytokinins (CKs) and strigolactones (SLs), which would relay the auxin signal by migrating into the bud. We developed a computational model of this system and designed experiments to identify and quantify sucrose action in this model. Experiments consisted in cultivating nodal stem segments of rosebush in vitro with different sucrose and auxin levels, and in quantifying bud elongation, CK level, and the expression of genes involve in SL biosynthesis and signalling.

Results. We observed that increasing sucrose level decreased the inhibition of bud elongation by auxin, so that buds fed with high sucrose level were less inhibited by a given amount of auxin than those fed with low sucrose level. In accordance with literature, auxin repressed CKs and stimulated the expression of SLs biosynthesis genes. We demonstrate that sucrose repressed SL signalling. From these results, we designed a model that reflects the combined action of sucrose and auxin on bud elongation. We validated it for its capacity to predict the effect of external CK supply for different sucrose levels.

Conclusions. We developed the first "bud" model accounting for the combined action of auxin and sucrose on bud burst. Initially observed for rosebush, these results were also validated in pea, demonstrating model genericity. Next step is to integrate the "bud" model at the scale of a plant to understand the role of sugars, together with hormones, in the spatio-temporal regulation of bud burst. In particular, we want to evaluate if sugar availability can explain the variations of bud burst observed in response to light environment.

Modelling physiological regulation of bud burst by sucrose and auxin

	Jessica Bertheloot, François Barbier, Frédéric Boudon, Maria-Dolores Perez-Garcia, Thomas
	Péron, Christine Beveridge, Christophe Godin, Soulaiman Sakr
	jessica.bertheloot@angers.inra.fr
	Research focus. Branching is an important process for productivity (number of productive
	branches) and for visual quality of ornamental plants (branches spatial arrangement).

Péron, et al.. Modelling physiological regulation of bud burst by sucrose and auxin. IEEE International Conference on Functional-Structural Plant Growth Modeling, Simulation, Visualization and Applications (FSPMA 2016), Nov 2016, QingDao, China. . hal-01400025


Salah Eddine Saidi
email: salah-eddine.saidi@ifpen.fr

Nicolas Pernet
email: nicolas.pernet@ifpen.fr

Yves Sorel
email: yves.sorel@inria.fr

Abir Ben Khaled
email: abir.ben-khaled@ifpen.fr

Acceleration of FMU Co-Simulation On Multi-core Architectures

Keywords: FMI, co-simulation, multi-core, scheduling, heuristic

The design of cyber-physical systems is a complex process and relies on the simulation of the system behavior before its deployment. Co-simulation allows system designers to simulate a whole system composed of a number of interconnected subsystems. Traditionally, these models are modeled by experts of different fields using different tools, and then integrated into one environment to perform simulation at the system-level. This results in complex and heavy co-simulations and requires adequate solutions and tools in order to reduce the execution time. Unfortunately, most modeling tools perform only mono-core simulations and do not take advantage of the omnipresent multi-core processors. This paper addresses the problem of efficient parallelization of cosimulations. It presents a multi-core scheduling heuristic for parallelizing FMI-compliant models on multi-core processors. The limitations of this heuristic are highlighted and two solutions for dealing with them are presented. The obtained speed-up using each of these solutions is illustrated and discussed for further improvements.

Introduction

Cyber-physical systems incorporate a combination of computational elements which collaborate in order to control physical processes. The complex nature of such systems requires cost, time and effort-effective design methodologies; therefore predicting their behavior and functioning scenarios before testing the real system is becoming more and more an indisputable step. Cosimulation aids in achieving these requirements as it allows the assessment of the design of the system by imitating its behavior. It consists mainly in simulating, on a computer, the global behavior of a multi-physics system composed of a number of interconnected subsystems. System designers can then identify potential design flaws and correct them before deploying the system.

Co-simulation faces however a number of challenges.

Actually, the simulated system is described by several interacting models which are often developed by experts of different fields using different tools and following different design approaches. The diversity of modeling tools and involved teams makes the coupling of the models a complex task. In fact, co-simulation necessitates efficient synchronized communications between the models where each model must be able to detect and respond to events of other models. Thanks to the FMI (Functional Mock-up Interface) standard [START_REF] Blochwitz | The Functional Mockup Interface for tool independent exchange of simulation models[END_REF], it is now possible to easily couple diverse models originating from different developers and tools. Nevertheless, executing FMI-compliant models raises some issues, which unless well handled, may reduce the cosimulation performance and limit the benefits of FMI. One major issue is the question of how to reduce the co-simulation execution time. Integrating heterogeneous models into one environment usually results in a complex and heavy to execute co-simulation which increases the demand of processing power.

As is well-known, increasing CPU frequency by means of silicon integration has reached its possible limits and semiconductor manufacturers switched in last years to building multi-core processors, i.e. integrating multiple processors into one chip allowing parallel processing on a single computer. Multi-core processors can reduce the execution time of a computational task by dividing it into several subtasks and assigning a subset of subtasks to each core to be processed in parallel. Most simulators, however, have mono-core simulation kernels and do not take advantage of the computation power brought by multi-core architectures. Therefore, enabling parallel execution of heavy co-simulations on multi-core processors is keenly sought by the developers and the users of simulation tools. However, fulfilling this objective is not trivial and appropriate parallelization schemes need to be applied on co-simulation models in order to accelerate their execution on multi-core processors.

FMI gives information about inputs and outputs relationships inside a model that is exported as a FMU (Functional Mock-up Unit). A model FMU is a package that encapsulates a XML file containing among other data the definitions of the model's variables, and a library defining the equations of the model as C functions. Input, output and state variables are updated by what we name "operations" which may call different functions provided by the FMU.

Given these features, various execution possibilities can be realized and the parallelization of co-simulation models on a multi-core processor can be seen as the following problem: Find an allocation of the different operations to the different cores and define an execution order, i.e. schedule the operations that are allocated to each core. When solving this problem, the utilization of the available cores has to be optimized in order to achieve the best acceleration. Using parallel computing terminology, the problem consists in finding a schedule for all the operations of the co-simulation on a multi-core processor. This paper deals with the problem of scheduling operations of heavy complex co-simulation models on multicore processors in order to accelerate the simulation execution. It follows the approach presented in [START_REF] Ben Khaled | Fast multi-core co-simulation of cyber-physical systems: Application to internal combustion engines[END_REF] by addressing two limitations of the previous work. First, an efficient multi-core scheduling can not be obtained without taking into account a good estimation of each operation's execution time. Second, the non-thread-safe implementation of FMUs prevent full exploitation of the potential parallelism of co-simulation graphs. Techniques for dealing with these limitations are here compared.

The rest of the paper is organized as follows. Next section presents related work on multi-core execution of simulations. Then our parallelization approach, firstly presented in [START_REF] Ben Khaled | Fast multi-core co-simulation of cyber-physical systems: Application to internal combustion engines[END_REF], is described in section 3, including a discussion about its present limitations. The fourth section presents our contribution, including the use of a toolchain for profiling cosimulation graph parallelism and explores the theoretical gain in execution speed-up over different architectures. Theoretical results are discussed and compared to real co-simulation executions in xMOD1 . xMOD is a co-simulation and a virtual experimentation platform, which allows mixing stand-alone and tool coupling cosimulations and the optimization of complex models execution. It provides a user-friendly interface in order to extend the simulation use to non-experts and ensure the continuity from Model-in-the-Loop to Hardware-in-the-Loop simulations. The last section concludes the paper and gives an outlook into our ongoing and future work.

Related Work

In order to achieve simulation acceleration using multicore execution, different approaches are possible and were already explored. From a user point of view, it is possible to modify the model design in order to prepare its multi-core execution, for example by using marked functions or Modelica extensions as in [START_REF] Elmqvist | Automatic gpu code generation of modelica functions[END_REF][START_REF] Gebremedhin | A data-parallel algorithmic modelica extension for efficient execution on multi-core platforms[END_REF]. From a modeling tool provider point of view, if providing OpenMP ready libraries is possible, the key feature for simulation acceleration is to provide techniques which offer speedup whatever the model is. Proposing parallel solvers or automatic parallel executions of model equations as in [START_REF] Elmqvist | Parallel model execution on many cores[END_REF][START_REF] Sjölund | Towards efficient distributed simulation in modelica using transmission line modeling[END_REF] is also an efficient way. In this paper, we address the problem from a co-simulation tool provider point of view. In such a tool, the user connects different FMUs, embedding solvers or not. In this case, it is not possible to change the models, the solvers, or the modeling tools. Such FMU assembly defines a graph of operations and the main opportunity to improve the co-simulation execution is consequently to accomplish an automatic parallelization of this graph. As shown in [START_REF] Ben Khaled | Multicore simulation of powertrains using weakly synchronized model partitioning[END_REF], splitting a model into several FMUs, by isolating discontinuities, may reduce the simulation time, even in the case of a mono-core execution. [START_REF] Ben Khaled | Fast multi-core co-simulation of cyber-physical systems: Application to internal combustion engines[END_REF] presented the RCOSIM approach. It consists in using each FMU information on input/output causality to build a graph, with an increased granularity and then exploiting the potential parallelism by using a heuristic to build an off-line multi-core schedule. This method has been tested on a real industrial model and significant speedup was obtained. This approach was implemented in the co-simulation tool xMOD and is available in its 2015 release.

3 Parallelization approach

Principle

The parallelization concept of xMOD is based on a task DAG (Directed Acyclic Graph) scheduling approach. Thanks to FMI, it is possible to access information about the internal structure of a model encapsulated in a FMU. In particular, FMI allows the identification of Direct Feedthrough and Non Direct Feedhrough outputs of a model. Since connections between different models of the co-simulation are also known, all data dependencies between the operations are known. Figure 1 shows an example of two models and their inter and intra-model dependencies.

The co-simulation can be described by a DAG where each vertex represents one operation and each edge describes a dataflow and thus a precedence constraint between two operations. The approach proceeds in two steps: First, the co-simulation DAG is constructed and then, the operations are allocated to the available cores in such a way to minimize the makespan of the graph. The makespan corresponds the execution time of the whole DAG.

The transformation of each model into an operation graph allows the parallelization of the model instead of considering it as an atomic block. Consequently the potential parallelism of the entire co-simulation is increased and can be better adapted to the hardware parallelism (number of cores in the case of a multi-core processor). The potential parallelism of a graph corresponds to vertices that are not dependent which characterize the partial order of the graph.

Model B Y A1 U B1 Y B2 U A1 Model A U A2 Y A2 Y A3 Y B1 U B3 U B2 Figure 1.
Inter and intra-model dependencies of two models.

Multi-core scheduling heuristic

The co-simulation DAG is built by exploring the relations between the models and between the operations of the same model. 2 illustrates the graph constructed from the two models of Figure 1.

In order to achieve fast execution of the co-simulation on a multi-core processor, an efficient allocation and scheduling of the DAG vertices has to be performed. xMOD uses an off-line scheduling heuristic similar to the one proposed in [START_REF] Grandpierre | Optimized rapid prototyping for real-time embedded heterogeneous multiprocessors[END_REF]. [START_REF] Ben Khaled | Fast multi-core co-simulation of cyber-physical systems: Application to internal combustion engines[END_REF] presented the use of this heuristic and the speed-up obtained by applying it on an industrial combustion engine model. The heuristic considers the execution time of each operation and aims at computing a schedule that minimizes the makespan of the graph.

Using the execution time C i of each operation OP i , the heuristic computes first the earliest start and end dates from the graph start denoted S i and E i , then the critical path CP := max E i (Algorithm 1). After that, the latest start and end dates from the graph end denoted S i * and E i * and then the flexibility At each step, the heuristic computes for a given operation the schedule pressure on a specific core. The sched-

F i = CP -E i -E i * are
Y A1 U B1 Y B2 U A1 U A2 Y A2 Y A3 Y B1 U B3 U B2 X B X A .
. ule pressure is the difference between the makespan increase, by allocating this operation to this core, and the operation's flexibility. The heuristic updates the set of candidate operations to be scheduled at each step. An operation is added to the set of candidate operations if it has no predecessor or if all of its predecessors have already been scheduled. The set of candidate operations holds the partial order associated to the graph. Then, for each candidate operation, the schedule pressure is computed on each core and the operation is allocated to its best core, the one that minimizes the pressure. After this step, a list of candidate operation-best core pairs is obtained. Finally, the operation with the largest pressure on its best core is selected and scheduled. The heuristic repeats this procedure until all operations are scheduled (Algorithm 3). This heuristic has originally been used to implement critical hard real-time applications where the execution times are usually estimated as the WCET (Worst Case Execution Time). On the contrary, co-simulation is not safety critical and the main goal here is to achieve fast execution, so average computation times can be used. Previously, execution times in xMOD were estimated based on the observation of practical examples as follows: update state operations are by far more costly so they are assigned significantly higher execution times then update out put operations, whereas update input operations are just data copy whose cost is negligible. The work presented in this paper attempts to improve the efficiency of the scheduling heuristic by measuring more accurately the execution times using a profiler as described in section 4.1.

Limitations of the approach

Although the presented scheduling heuristic resulted in interesting co-simulation speed-ups, it has some limitations that have to be considered in the multi-core scheduling problem in order to obtain better performances. First, so far, the multiprocessor scheduling heuristic uses empiric operations execution times. By using realistic execution times for each operation, the multi-core execution of the simulation should be improved. In this paper, we present some results, based on a profiling technique.

Second, FMI standard does not presently require that FMU functions have to be thread-safe, i.e. they cannot be executed simultaneously as they may share some resource (variables for example) that might be corrupted if two operations try to use it at the same time. This implies that at any instant during the execution of the co-simulation, one and only one operation of the same FMU can be executed. Consequently, if the scheduling heuristic allocates two or more operations belonging to the same FMU to different cores, a mechanism that ensures these operations are executed in strictly different time intervals must be set up.

This section presents a theoretical study of the achievable speed-up on a use-case, using the SynDEx 2 software [START_REF] Sorel | Syndex: System-level cad software for optimizing distributed real-time embedded systems[END_REF][START_REF] Sorel | From modeling/simulation with scilab/scicos to optimized distributed embedded real-time implementation with syndex[END_REF]. Then, these theoretical results are compared with xMOD co-simulation runs, with two different implementations for guaranteeing a mutual exclusion between different operations of the same FMU.

Toolchain

A toolchain is proposed to assist the developer in parallelizing co-simulations. Using this toolchain, it is possible to assess new solutions before implementing them in xMOD thanks to the SynDEx software. SynDEx is a system level CAD software based on the Algorithm-Architecture Adequation (AAA) methodology [START_REF] Sorel | Real-time embedded image processing applications using the algorithm architecture adequation methodology[END_REF]. It was developed to optimize the implementation of real-time distributed applications onto multicomponent architectures. The workflow is illustrated in Figure 3. When different FMUs are imported into xMOD and connected together, a file which describes inter-model connections is generated. This file and the XML files of the different FMUs of the co-simulation are passed to a converter which parses the files and produces equivalent files (.sdx) compliant to the SynDEx format. The co-simulation code is profiled in order to obtain the execution times of the different operations which are introduced in SynDEx. SynDEx offers the possibility to use the multi-core scheduling heuristic outlined in this paper , as well as other kinds of heuristics, and therefore makes it possible to study the achievable co-simulation speed-up before implementing the heuristic in xMOD.

Use-case description

In this work, experiments have been carried out on a Spark Ignition (SI) RENAULT F4RT engine cosimulation using 5 FMUs. It is a four-cylinder in line Port Fuel Injector (PFI) engine in which the engine displacement is 2000 cm 3 . The air path is composed of a turbocharger with a mono-scroll turbine controlled by a waste-gate, an intake throttle and a downstreamcompressor heat exchanger (Figure 4). The engine model was developed using ModEngine library [START_REF] Benjelloun-Touimi | From physical modeling to real-time simulation: Feedback on the use of Modelica in the engine control development toolchain[END_REF]. ModEngine is a Modelica library that allows for the modeling of a complete engine with diesel and gasoline combustion models. The engine model was imported into xMOD using the FMI export features of the Dymola tool. This use-case has over 100 operations which are scheduled by the multi-core scheduling heuristic.

Results and Discussions

Using the toolchain, a .sdx file of the use-case was generated in order to evaluate, in SynDEx, the theoretical speed-up obtained by parallelizing the model on different numbers of cores, using the multi-core scheduling heuristic of section 3.2. For each schedule the speedup is computed by dividing the mono-core schedule makespan by the schedule makespan. Figure 5 gives the different theoretical speed-ups. The best speed-up is close to 3,6 and is reached with 6 cores. Finding the minimal number of cores which offers the maximum speed-up is interesting if a large number of simulation runs (possibly with different parameters) have to be performed: If a large number of cores is available, multiple runs could be launched in parallel with the adequate number of cores dedicated to each run. This research of the minimal number of necessary cores to reach the maximum speed-up could be scripted and automatically performed before the simulation.

In order to tackle the constraint of non thread-safe FMU functions, two mutual exclusion strategies have been implemented in xMOD and the performance obtained using each of them has been evaluated. The first one does not modify the multi-core scheduling heuristic result and uses a dedicated mutex (system object that guarantees mutual exclusion) for each FMU: Every time a FMU function call is made at runtime, the associated mutex have to be acquired before the execution of the function code can be started. The second solution is explained in [START_REF] Ben Khaled | Fast multi-core co-simulation of cyber-physical systems: Application to internal combustion engines[END_REF] and consists in modifying the multi-core scheduling heuristic to always allocate the operations of a same FMU to the same core (constrained allocation). If constrained allocation is used, the search space of the scheduling heuristic is reduced, i.e. at each step, for a given candidate operation, if there is another operation of the same FMU that has already been allocated to a specific core, the candidate operation is allocated to this same core without the need to test it on the other cores. Thanks to SynDEx, it is easily possible to theoretically estimate the impact of using the constrained allocation in the multi-core scheduling heuristic. Results are given in Figure 5. It shows that the expected speed-up in the case of constrained allocation is less than the one using unconstrained allocation, when the number of cores is less than 5, but similar when 5 cores or more are available. When using less than 5 cores, the large number of update out put operations can be efficiently allocated only if the unconstrained allocation is used: The speed-up difference between the constrained and unconstrained allocation cases is due to this restriction on the allocation. Five is the minimal number of cores for enabling the execution of each update state operation on a different core. Due to the predominant execution times of the update state operations, their impact on the speedup overrides the possibility of optimizing the allocation of the other operations. This explains why the speedup difference between the unconstrained and the constrained allocation cases becomes very small with 5 cores or more. In order to compare the two mutual exclusion strategies, we implemented them in xMOD. Execution times measurements were performed by getting the system time stamp at the beginning of the simulation and after 30 seconds of the simulated time. As previously, we compare the speed-up by dividing the mono-core simulation execution time by the simulation execution time on a fixed number of cores. Figure 6 sums up the results, where unconstrained allocation corresponds to the use of mutex objects. It shows the impact of mutex overhead on the speed-up. Whatever the number of the available cores, the speed-up remains close to 1,3. On the contrary, the implementation in xMOD of the constrained allocation gives similar results in terms of speedup improvement when increasing the number of cores until 5. Nevertheless, the maximum measured speed-up (2,4) remains smaller than the theoretical one (3,5). In fact, the theoretical speed-up computation considers the makespan ratio without estimating any synchronization cost between cores. The real implementation in xMOD contains synchronization objects between operations to ensure the consistency of data dependencies which certainly have an important impact on the speed-up.

Conclusion and Future Work

The work presented in this paper dealt with the problem of co-simulations acceleration by means of parallelization on multi-core processors. We proposed to extend our previous work by taking into account operations execution times in the multi-core scheduling heuristic. This allows the optimization of the number of the dedicated cores to the simulation, by performing architecture exploration with SynDEx. Our experiments in xMOD on an industrial use-case, gave important speed-up results (2,4). Nevertheless, it also shows the impact of the mutual exclusion constraint on the co-simulation acceleration. Providing thread-safe FMU implementation could offer important simulation acceleration opportunities. In our ongoing work, we are exploring graph transformation techniques to improve the handling of the mutual exclusion problem of FMUs. We also envision to extend these results to the multi-rate co-simulation of FMUs by developing an efficient multi-core scheduling heuristic to handle it.

Figure 2 .

 2 Figure 2. Dependency graph of the models of Figure 1.

  Figure 4. Spark Ignition (SI) RENAULT F4RT engine model.

Figure 3 .Figure 5 .

 35 Figure 3. Proposed toolchain to assist in the development and assessment of scheduling heuristics.

Figure 6 .

 6 Figure 6. Measured speed-up.

  out put operations and the update state operation of the same model, because the computation at instant k of an output Y k may need the state variable X k which is no longer available after update state computed X k+1 . Running the cosimulation consists in executing the graph repeatedly. At each co-simulation step the whole graph is executed and a new execution of the graph cannot be started unless the previous one was totally finished because some outputs at the current step may depend on certain state variables computation at the previous step. Such outputs are not identified in the model's FMU and thus all outputs are considered as dependent on the state variables. Figure

The operations are either update out put , update input or update state . A vertex is created for each operation and edges are then added between vertices if a data dependency exists between the corresponding operations. This information can be extracted from the model's FMU. Since FMI does not give information about the dependencies between the state variables computation and the input and output variables computations, edges connect all update input operations and the update state operation of the same model, which means that all inputs are needed in order to update the state. Furthermore, edges are placed between all update

  computed (Algorithm 2). Computation of S i , E i and CP

	Initialization;
	Set Ω the set of all the operations;
	Set O the set of operations without predecessors;
	foreach OP i ∈ O do
	S i := 0; E i := S i +C i ;
	end
	Set O the set of operations whose all immediate predecessors
	were treated;
	while O = / 0 do
	foreach OP i ∈ O do
	S end
	end
	CP := 0;
	foreach OP i ∈ Ω do
	if CP < E i then
	CP := E i ;
	end
	end
	Algorithm 1: Initialization;
	Set Ω the set of all the operations;
	Set O the set of operations without successors;
	foreach OP i ∈ O do E * i := 0; S * i := E * i +C i ;
	end
	Set O the set of operations whose all immediate successors were
	treated;
	while O = / 0 do
	foreach OP i ∈ O do
	E * i := max(S * h : OP i → OP h );
	(OP h are the immediate successors of OP i ); S * i := E * i +C i ; Remove OP i from the set O ;
	Add to the set O all predecessors of OP i for which all
	successors were already scheduled;
	end
	end
	foreach OP i ∈ Ω do
	F

i := max(E h : OP h → OP i ); (OP h are the immediate predecessors of OP i ); E i := S i +C i ; Remove OP i from the set O ; Add to the set O all successors of OP i for which all predecessors were already scheduled; i := CP -E i -E * i ; end Algorithm 2: Computation of S * i , E * i and F i

  Core j ); (new start date of OP i when executed on Core j ); cost i, j := S i +C i + E * i -CP; (cost of OP i when executed on Core j ); if cost i, j < cost i then Set cost i := cost i, j ; Set BestCore i := Core j ; end end end Find OP i with maximal cost i in O; Schedule OP i on its core BestCore i ; Set k := BestCore i ; T Core k := T Core k +C i ; (Advance the time of Core k ); Remove OP i from the set O; Add to the set O all successors of OP i for which all predecessors are already scheduled;

	Initialization;
	Set Ω the set of all the operations;
	Set Γ the set of all the available cores;
	Set O the set of operations without predecessors;
	while O = / 0 do
	foreach OP

i ∈ O do Set cost i to ∞; (cost of OP i is set to the maximum value); foreach Core j ∈ Γ do S i := max(S i , T end

Algorithm 3: Multi-core scheduling heuristic

www.xMOD-software.com