Roman Iakymchuk
email: roman.iakymchuk@telecom-sudparis.eu

François Trahay
email: francois.trahay@telecom-sudparis.eu

LiTL: Lightweight Trace Library

Keywords: Performance analysis, binary trace library, multithreaded applications, hybrid programming models

In order to design parallel applications that exploit efficiently modern supercomputers with their complex architecture, it is required to apply a variety of programming models and techniques. The development of such applications as well as their tuning to achieve the optimal performance is practically impossible without performance analysis tools.

In this paper, we present a scalable binary trace library -Lightweight Trace Library (LiTL). LiTL aims at providing performance analysis tools with a low-overhead event recording service. In order to enhance the scalability and the performance of such tools, we implement various optimization strategies and techniques in LiTL. The experiments on a suite of synthetic and standard computation kernels show that the overhead caused by LiTL in conjunction with EZTrace remains low. Furthermore, LiTL often performs better than the existing event recording libraries, saving up to 88 % of the CPU time and up to 68 % of the disk space. Finally, LiTL is a generic library that can be used in conjunction with a wide range of performance analysis tools.

I. INTRODUCTION

High-performance computing (HPC) is no longer a topic related only to Computer Science. Scientists from many other disciplines such as Physics, Geodesy, and Biology rely more often on parallel computing for their experiments. In general, HPC clusters permit to run larger simulations and/or obtain results in a rather short time. However, due to the complexity in hardware -multi-core processors, NUMA architectures, GPUs accelerators -as well as the hybrid programming models that mix MPI, OpenMP, and/or CUDA, the development of a scalable application that exploits efficiently HPC platforms is a difficult task.

Performance analysis tools provide a great help during the development, tuning, and optimization of parallel applications. The purpose of such tools is to depict the application behavior in order to find performance bottlenecks and improve the applications performance. Those tools are usually used to gather the information about the program execution for immediate or further performance analysis. One way of gathering such information is to store time-stamped events in a trace file. Those events represent the important phases of the application execution such as calls to MPI communication primitives and synchronization among threads. Collecting this data requires to intercept the calls to the corresponding functions and to record events in a trace file.

In order to depict accurately the behavior of an application, it is crucial that the performance analysis tool does not disturb the program execution: the overhead of such tool should be as small as possible. Previously, we presented a lightweight mechanism for instrumenting functions [START_REF] Aulagnon | Runtime function instrumentation with EZTrace[END_REF]. Once functions have been intercepted, the performance analysis tool needs to record events in a trace file with the lowest overhead possible. In this paper, we focus on the second important component of performance analysis tools: the event recording mechanism.

We present LiTL -a lightweight binary trace library. LiTL aims at providing performance analysis tools with a scalable event recording mechanism that utilizes minimum resources of the CPU and memory. In order to efficiently analyze modern HPC applications that combine OpenMP (or Pthread) threads and MPI processes, we design and implement various mechanisms to ensure the scalability of LiTL for a large number of both threads and processes. Although LiTL is a generic library that can be used in conjunction with many performance analysis tools, here we integrate it in the EZTrace framework as a replacement for FxT [START_REF] Danjean | An efficient multi-level trace toolkit for multi-threaded applications[END_REF] -the current event recording library. The evaluation shows that LiTL performs better than FxT for multi-threaded applications. LiTL permits to analyze larger MPI jobs. Moreover, trace files that are generated by LiTL are smaller in size than the ones created by FxT.

The rest of the paper is organized as follows. Section II provides literature overview of existing performance analysis tools as well as trace libraries. Section III presents the EZTrace framework and its core components including the current event recording library. The main features of the new scalable binary trace library are introduced in Sections IV to VI. We evaluate LiTL in conjunction with EZTrace and compare it with FxT in Section VII. Finally, Section VIII discusses conclusions and the future work.

II. LITERATURE OVERVIEW

When developing high-performance or parallel applications, a big concern is put on the communication between threads/processes/nodes and also the performance optimization. In order to create scalable software with the optimal execution time, software developer often rely on performance analysis tools. Numerous research has been conducted on this topic. As a result, many tools were designed for tracing the execution of parallel applications to better understand their behavior as well as to find performance bottlenecks. Some tools are dedicated to a specific programming model or library, such as MPI Parallel Environment (MPE) [START_REF] Chan | An efficient format for nearly constant-time access to arbitrary time intervals in large trace files[END_REF] that targets MPI applications; POSIX Thread Trace Toolkit (PTT) [START_REF] Decugis | NPTL Stabilization Project[END_REF] aims at Pthread-based applications; OMPtrace [START_REF] Caubet | A Dynamic Tracing Mechanism for Performance Analysis of OpenMP Applications[END_REF] instruments OpenMP applications. Other performance analysis tools cover multiple models and can track calls to multiple libraries on various architectures, for instance VampirTrace [START_REF] Müller | Developing scalable applications with vampir, vampirserver and vampirtrace[END_REF], TAU [START_REF] Shende | The Tau Parallel Performance System[END_REF], Paraver [START_REF] Pillet | Paraver: A tool to visualize and analyze parallel code[END_REF] or Scalasca [START_REF] Geimer | The Scalasca performance toolset architecture[END_REF]. These tools can be used for instrumenting custom libraries or applications through manual or automatic instrumentation of their code.

The performance analysis tools that generate execution traces rely on tracing libraries dedicated to a specific trace format. Some formats are limited to a particular programming model. For instance, the Open Trace Format (OTF) [START_REF] Knüpfer | Introducing the open trace format (otf)[END_REF] focuses on MPI events and function entry/exit. OTF2 [START_REF] Eschweiler | Open trace format 2-the next generation of scalable trace formats and support libraries[END_REF] is an extension of OTF that adds support for OpenMP events. Some other trace formats are self-defining, allowing to depict any event in a more generic way. For instance, the Pajé file format [START_REF] De Oliveira Stein | Pajé trace file format[END_REF] and the Pablo self-defining data format. The concept of the other trace formats is to store data in a raw format within a trace file. The application that reads the trace file is then in charge of adding meaning to the recorded events. The FxT [START_REF] Danjean | An efficient multi-level trace toolkit for multi-threaded applications[END_REF] library create traces using the raw format.

Additionally to the optimization applied in file formats, some libraries -such as OTF or OTF2 -use a compression mechanism to reduce the file size. This is done by using a general purpose compression library like zlib. The main drawback of file compression is that it causes an additional overhead when writing the traces to the disk. We do not consider the trace compression here, however adding the compression phase before writing traces to the disk is rather straightforward.

III. EZTRACE

EZTrace is a general framework for analyzing highperformance and parallel applications. EZTrace relies on a tracing mechanism that aims at executing an application once in order to record its execution trace. After the execution of the application, EZTrace analyzes the resulting trace and extracts various statistics. Furthermore, EZTrace provides a functionality to generate trace files in different formats. These files can be visualized with tools such as Vampir [START_REF] Müller | Developing scalable applications with vampir, vampirserver and vampirtrace[END_REF], ViTE [START_REF] Coulomb | ViTE -Visual Trace Explorer[END_REF], or Triva [START_REF] Schnorr | Triva: Interactive 3d visualization for performance analysis of parallel applications[END_REF].

EZTraceis shipped with modules for the main libraries that are used in parallel programming, e.g. MPI and OpenMP, and it allows third-party developers to create modules that instrument their own functions [START_REF] Trahay | EZTrace: a generic framework for performance analysis[END_REF]. For instance, the PLASMA [START_REF] Agullo | Numerical linear algebra on emerging architectures: The PLASMA and MAGMA projects[END_REF] linear algebra solver integrates an EZTrace module that permits to visualize the tasks scheduling operated by PLASMA.

For recording trace files, EZTrace modifies the binary of the application on-the-fly and inserts probes before and after functions of interest [START_REF] Aulagnon | Runtime function instrumentation with EZTrace[END_REF]. Then, the program is executed in the regular way. When the program enters one of the instrumented functions, e.g. MPI_Send, EZTrace records an event (MPI_SEND_ENTRY), executes the actual function (MPI_Send), and records another event that marks the exit of the function (MPI_SEND_EXIT). Each event consists of a timestamp, an event code, a thread identifier, and additional data depending on the event type. In case of the MPI_SEND_ENTRY event, EZTrace records the destination process, the message length, and the message tag as the additional data. At the end of the application execution, events are written into trace files.

After the recording phase, EZTrace reads the trace files that were generated and analyzes the events. Depending on the user choice, during the analysis phase various statistics, e.g. the time spend in the MPI communication, can be extracted and Pajé or OTF files can be generated in order to visualize the application execution.

A. Recording events with FxT

For recording events, EZTrace relies on the FxT library. FxT is a multi-level trace library [START_REF] Danjean | An efficient multi-level trace toolkit for multi-threaded applications[END_REF] that permits to record events in a binary format from both kernel-and user-space. Even though FxT is efficient for small-scale applications, it suffers from several scalability issues.

The first scalability issue is related to the event recording mechanism: FxT generates one trace file per process during the execution of the application. Moreover, the recording mechanism fails to scale when the number of threads per process grows. For instance, on multi-core architectures each process may contain dozens of threads that write events into a single event buffer -increasing both the contention and the application execution time. Thus, the FxT scalability issue concerning the large number of threads is a big obstacle on modern CPUs.

Afterwards, during the post-mortem analysis all the recorded trace files are read by EZTrace for further analysis. Due to the constrains imposed by the Linux operating system, there is a limitation on the number of trace files (it typically equals 1024 files) that can be opened at once by the process, which is in charge of the performance analysis. Thus, EZTrace currently cannot be used for large MPI jobs.

Performance analysis may impact heavily the efficiency of parallel applications, especially the large scale ones that generate millions of events producing huge trace files. A major issue when analyzing such applications is the size of generated trace files. Reducing the size of traces is critical in order to make the analysis of large scale applications lightweight and harmless for their performance. Although compression techniques can be used, it is also crucial to reduce the size of data to its maximum before compression.

Finally, the complex multi-level tracing mechanism that permits to record events from both kernel-and user-space makes FxT difficult to maintain. Since EZTrace records events only from user-space, a large portion of the source code as well as the data structures are redundant. Due to the complex structures and the functionality behind, adding new features that would improve the scalability of FxT is an arduous task that is comparable to the development of a new trace library with the same functionalities from scratch.

B. Recording events with LiTL

In order to resolve the above-mentioned issues, we propose to design a new trace binary library that provides similar functionalities to FxT:

• Records events from user-space with a low overhead;

• Improves the scalability, especially for applications that rely on a large number of threads or processes.

The new library is named Lighweight Trace Library (LiTL). LiTL can be used by performance analysis tools such as EZTrace for recording events in trace files. As in FxT, LiTL minimizes the usage of the CPU time and memory space in order to avoid disturbing the application that is being analyzed.

LiTL is also fully thread-safe that allows to record events from multi-threaded applications. In contrast to FxT, LiTL records events only from user-space. This permits to apply only the relevant data structures and functionalities, and thus to simplify the maintainability of the library.

In the next Sections, we present various techniques that permit to resolve the scalability issues of FxT. Section IV describes the LiTL approach for reducing the size of the generated trace files. Section V depicts the strategy for improving the scalability of the library when a large number of threads concurrently record events. Finally, Section VI presents the solution that allows LiTL to process large MPI jobs. All of these features are already implemented in LiTL.

IV. UTILIZED RESOURCES

EZTrace gathers the information regarding the program execution in form of events that are generated before and after calls to the instrumented routines.

As an example, let us consider a simple MPI ping-pong program that relies upon MPI_SEND and MPI_RECV routines. These routines have the following prototypes: We develop an MPI module for EZTrace, where we determine a list of key functions that are instrumented. For each function, we also define a list of parameters that are recorded and can be useful during the post-mortem analysis. For instance, for MPI_SEND we record two events as follow:

1 EZTRACE_EVENT4(EZTRACE_MPI_START_SEND, size, dest, tag, comm); 2 EZTRACE_EVENT3(EZTRACE_MPI_STOP_SEND, dest, tag, comm); EZTRACE_MPI_START_SEND is recorded before the call to MPI_SEND and EZTRACE_MPI_STOP_SEND afterwards. Also, in these events we store the information regarding four parameters, namely a message size, a message tag, a destination, and an MPI communicator. This procedure is applied for all intercepted functions.

Each event in the FxT library consists of two parts: the event core (the event code, the time when the event occurred, the thread identifier, and the number of parameters) and event parameters. The number of event parameters recorded by FxT varies from zero to six. Also, there is an additional event type in FxT, where parameters are stored in the string format. Its purpose is to gather either the regular parameters in a string format or the information about the abnormal behavior of applications like thrown exceptions.

As it can be seen from the example above, the parameters passed to each event have different data type. In order to handle the variety of possible cases, event's parameters in FxT are represented by the largest data type, which is uint64_t on x86_64 architectures. Hence, any parameter -no matter whether it is a char, an int or a long int -can be recorded without being truncated. However, the reserved slot for each parameter is often bigger than its actual size. Thus, this leads to the non-optimal usage of resources. For instance, in case of MPI_SEND each parameter requires 4 bytes of the storage, but it occupies 8 bytes instead.

From the experiments we observe that in many cases the required space to store an event is smaller than the one reserved by FxT. Hence, we apply our optimization strategy to utilize every byte from both the memory and disk space. Such optimization is crucial, especially for applications that generate millions of events and produce large trace files. Moreover, writing a huge amount of data to the disk increases the overhead of EZTrace. So, our goal is to keep trace files as small as possible without losing any of the recorded data. Therefore, we propose to use the compacted event storage that aims at utilizing every byte from the allocated space.

In our approach, we introduce three different types of events: regular, raw, and packed. The regular event is similar to the one in FxT, meaning without any major optimization being involved. The packed and raw (parameters are stored in a string format) events represent the optimized versions of storing events, where each parameter can be saved as a group of bytes. Accordingly, by using the event type packed for recording and storing events, we theoretically are capable to save up to 65 % of the disk space1 compare to both the regular LiTL and FxT events. Also, we increase the number of parameters to 10, which covers a wider range of possible cases.

Fig. 1(a) shows, on an example of three events with different number of parameters, the occupied space of events within the trace file recorded by EZTrace with FxT. We symbolically partitioned the trace file into bytes and also chunks of bytes, which store event's components. The space occupied by each event with its parameters is highlighted with parentheses. Fig. 1(b) shows the storage of the recorded events in the trace file while using EZTrace with LiTL. We consider one particular scenario when each event's parameter can be represented by uint8_t; this requires only one byte for the storage. To store larger event's parameters we use arrays of uint8_t. This scenario corresponds to the optimal performance in terms of the memory and disk space usage. Under this approach, not only the size of the core event's components is shrunk, but also the size of event's parameters is reduced significantly. The gained performance, e.i. the reduced space, can be characterized by the gray area that corresponds to the difference in storage between FxT and LiTL events. The size of three events recorded by LiTL is smaller than the size of one event with five parameters recorded by FxT. This figure confirms our assumption regarding the possibility of reducing the size of both the recorded events and trace files.

V. SCALABILITY VS. NUMBER OF THREADS

The advent of multi-core processor have led to the increase in the number of processing units per machine. It becomes usual to equip a typical high performance computing platform with 8, 16, or even more cores per node. In order to exploit efficiently such facilities, developers can use hybrid programming models that mix OpenMP (or Pthread) threads and MPI processes within one application. Hence, the number of threads per node, which executes the same application, can be quite large -8, 16, or even more threads. The number of threads per node is the scalability issue for FxT, because in its implementation all threads within one process record events into a single buffer, see Fig. 2(a).

The FxT recording mechanism causes a contention problem -when multiple threads record events simultaneously, the pointer to the next available slot in the buffer is changed concurrently. The modifications of the pointer can be done atomically in order to preserve the data consistency. However, the atomic operation does not scale quite well when it is performed by a large number of threads at the same time. Thus, analyzing OpenMP applications that run lots of threads using EZTrace with FxT may result in the high overhead.

This scalability issue can be clearly seen on NUMA architectures where the buffer in which threads record events is allocated on one memory bank of the NUMA machine. When multiple threads gather application's events from different NUMA nodes, the overhead of the event recording may increases due to NUMA effects. Therefore, analyzing the execution of applications on NUMA machines with a decent amount of cores is potentially harmful for their performance.

A. Recording Events

While designing LiTL, we aim at resolving the abovementioned limitation of FxT. Thus, we propose to record events into separate buffers, meaning to have one buffer per thread instead of one buffer per process. This approach is illustrated on Fig. 2(b).

To keep multiple buffers in order within the trace file, we add a header into the trace file with the information regarding the number of threads and pairs <tid, offset>; tid stands for the thread identifier; offset corresponds to the position of the first chunk of events for a given thread within the trace starting from its beginning. The last event of each chunk contains either an offset to the next chunk of events or a symbol specifying the end of recording for a given thread. While flushing the current buffer to the trace file, the following two actions are performed:

1) Setting the offset of the current chunk to specify the end of the recording; 2) Update the offset from the previous chunk to point to the current one.

Fig. 3 demonstrates the storage mechanism on an example of three threads, including the positioning of chunks of events as well as the way of linking those chunks into one chain of the corresponding thread using offsets.

During the application execution, it may occur that some threads start recording events later than others. This scenario requires appropriate modifications and adjustments to the above approach. According to the previous approach, the header is the first block of data that is added to the trace file; it is written before flushing the first chunk of events. Thus, the header contains the information only regarding the started threads. In order to add pairs <tid, offset> of the late threads, we reserve a space for 64 pairs (chunk of pairs) between chunks of events within the trace file. So, when one among those late threads wants to flush its buffer to the trace file, we add its pair <tid, offset> directly to the next free spot in the chunk of pairs. The chunks of pairs are binded with offset in the same way as chunks of events. Therefore, EZTrace does not have limitations on the number of threads per process and also processes.

B. Post-Mortem Analysis

We develop the functionality for analyzing the generated traces by capturing the procedure of the event recording H E A D E R LiTL_v OS #threads buf_size tid0 offset0 tid1 offset1 tid2 offset2 chunk00 events offset01 chunk20 events offset21 chunk01 events offset02 chunk10 events offset11 chunk02 events offset03 Fig. 3: Storage of events recorded by EZTrace with LiTL on multi-threaded applications. In the figure, LiTL_v contains information about LiTL; OS -about OS and architecture; #threads stands for the number of threads; buf_size -the buffer size. mechanism. At first, LiTL reads the trace header with the information regarding the buffer size, threads (the number of threads, tids, and offsets), and also pairs <tid, offset> that correspond to the late threads. Using this preliminary information, LiTL allocates memory buffers for reading; the number of buffers equals the number of threads used during the recording phase, meaning one buffer per thread. Then, LiTL loads chunks of events from the trace file into these buffers using pairs <tid, offset>. After processing the first chunks of events, LiTL loads the buffers with the next ones using the information concerning their positions in the trace, which is given by the offsets. This procedure is recursive and stops when the symbol specifying the end of recording is reached.

VI. SCALABILITY VS. NUMBER OF TRACES

In FxT, one trace file is generated per process. This means that for parallel applications with hundreds of MPI processes the equal amount of trace files is created. This is one side of the problem. The other side appears while analyzing the applications execution due to the limitation on the number of trace files that can be opened and processed at the same time. Therefore, FxT does not perform well and even crashes when the number of traces exceeds the Linux OS limit on the number of simultaneously opened files.

In order to overcome the opened files limitation imposed by the Linux OS, one may increase the limit to the maximum possible value. However, this would temporarily solve the problem. Instead, we propose to create archives of traces during the post-mortem phase. More precisely, we suggest to merge multiple traces into a trace archive using the litl_merge utility from LiTL. Fig. 4 illustrates the structure of the new combined trace created by litl_merge. The archives of traces preserve all information concerning each trace: headers, pairs <tid, offset>, and positioning of events chunks. They also contain new global headers that store the information regarding the amount of trace files in the archive and triples <fid, size, offset>; fid stands for a file identifier; size is a size of a particular trace file; offset holds the position of a trace file within the archive of traces. Therefore, archives of traces not only solve the FxT problem, but also make the further analysis of the applications performance more convenient. The archives can be read by litl_read and, therefore, analyzed by EZTrace.

One more useful feature provided by LiTL, which is the opposite of litl_merge, is a possibility to extract trace files from archives with the litl_split utility. This utility can be applied when there is a need to analyze a particular trace or a set of traces among the merged ones.

VII. PERFORMANCE RESULTS

When analyzing the execution of high-performance or parallel applications that generate millions of events, it is important to keep the overhead of the event recording library as low as possible. In this Section, we evaluate the performance and the overhead of both FxT and LiTL in conjunction with EZTrace on synthetic benchmarks as well as on standard application kernels.

A. Experimental Setup

We use the RWTH cluster that is composed of 60 nodes connected through DDR InfiniBand [START_REF] Mey | The RWTH HPC-Cluster User's Guide: Version 8.2.4[END_REF]. Each node is equipped with two quad-core Intel Xeon E5450 (Harpertown) processors and 16 GB of RAM. Each of the four cores on Harpertown operates at 3.0 GHz.

The RWTH cluster's nodes execute the Scientific Linux operating system release 6.1 (Carbon) with the 2.6.32 Linux kernel. There, we installed EZTrace v.0.9 with both FxT and LiTL, the NAS Parallel Benchmarks (NPB) v.3.3, and Multi-Zone version of NPB v.3.3 (NPB-MZ). These libraries were compiled with the gcc v.4.6.1 compiler under the optimization level three enabled. To run both NPB and NPB-MZ, we use MPI from the Open MPI v.1.6.4 library.

B. Utilized Resources

To validate our optimization strategy concerning the event size and the occupied disk space, we created a synthetic benchmark and conducted a set of tests recording

• Events with various number of parameters;

• Events with a homogeneous type of parameters in order to have a better picture regarding the size of each event at the end.

Each test case (an event with fixed number of parameters of one data type) were executed 10 6 times, meaning that 10 6 events were recorded.

The results of the conducted experiments are reported in Tab. I. FxT does not make a difference between various data types and reserves the maximum amount of space for each parameter. For instance, the FxT traces composed of events with six parameters of char, int, or double data types are equal in terms of size. In contrast, LiTL uses compact storage and its traces have the following characteristics

• Only the required space for the storage is used. For example, the difference in size between traces composed of events without parameters and the ones with one parameter of type char is exactly 10 6 . This corresponds to the size of 10 6 chars;

• The occupied storage varies depending on the data type. For instance, events with one character and one integer require 18 and 21 bytes of the storage, accordingly.

C. Scalability vs. Number of Threads

In this subsection, we verify our solution on improving the scalability of EZTrace and reducing the contention while analyzing the performance of multi-threaded applications. In FxT, multiple threads within a process record events into a single buffer. This causes contention and the performance overhead, in particular on NUMA architectures. Our solution for this problem is to use multiple buffers, more precisely one buffer per thread. To validate this solution, we develop a synthetic multi-threaded program that calls repeatedly a simple function. This function is intercepted by EZTrace and two events are recorded: the first is at the entry of the function and the second at its termination. Thus, when the function is called for numerous times, the program generates a large number of events. Since the computation are simple and fast, and the function is called by multiple threads, this test is a good model of event-intensive multi-threaded applications.

We evaluate the performance of EZTrace with both FxT and LiTL running the multi-threaded test program on the Harpertown nodes. For comparison, we also evaluate the performance of VampirTrace on the same program. Tab. II presents the average time that is required to record one event on the Harpertown nodes. When the number of threads used at the program execution increases, the time for recording an event with FxT growth rapidly. This is due to the contention between the threads that try to record events concurrently. With LiTL and VampirTrace, the average time increases, but not as dramatically as with FxT. FxT and LiTL deliver nearly the same performance for single-threaded applications. Nevertheless, the scalability of LiTL when the number of concurrent threads grows is much better, allowing LiTL to perform up to 8.4 times faster.

D. Scalability vs. Number of Traces

In order to verify the scalability of both event recording libraries with the number of recorded traces, we created a simple MPI program that generates few MPI events during its execution. We ran the program using 256 MPI processes on the RWTH cluster. As a result, 256 trace files were created. Afterwards, we extracted the statistic information (the size of messages, the time spent in the communication, and the number of events handled) from the generated traces using the eztrace_stats utility from EZTrace. At first, eztrace_stats performed well on the 256 traces for both FxT and LiTL. Then, we decreased the limit on the number of simultaneously opened files from 1024 to 200 using the ulimit command and repeated the experiment. At this time, eztrace_stats was not able to gather the statistic information with both event recording libraries. In order to overcome this limit, we applied litl_merge and combined traces into two archives; each of them stores equal amount of traces. Consequently, we succeeded with the performance analysis of the program using EZTrace with LiTL. The same technique is applicable for eztrace_convert, which produces traces in different formats. Moreover, EZTrace supports analysis of both traces and archives of traces mixed together.

E. NAS Parallel Benchmarks (NPB)

To evaluate the overhead of both FxT and LiTL within EZTrace on more realistic scenario, we measure their performance for NAS Parallel Benchmarks (NPB). NPB is a set of computation kernels derived from computational fluid dynamics applications to evaluate performance of parallel and high-performance applications. Later, the benchmark suite was extended including new kernels, e.g. for parallel I/O and computational grids. The problem sizes of those kernels are predefined and divided into eight classes; this value depends upon the kernel [START_REF]NAS Parallel Benchmarks[END_REF]. We conducted our experiments using 64 computing processes for problems from the class C, which are already decent to make the CPUs busy with the computations and also to generate a good amount of events. The MPI functions of the NPB kernels were instrumented with EZTrace and recording events was performed with both FxT and LiTL. Each benchmark was executed at 10 times and then the median time was selected. Tab. III presents the time measurements: the median execution time and the overhead of both FxT and LiTL from the raw NPB execution. The raw as well as the instrumented kernels show the fluctuation in the measurements. This is an explanation for the negative time overhead. Nevertheless, the timings show that instrumenting the NPB kernels with EZTrace using both FxT and LiTL has small or even negligible impact on the performance: for all experiments EZTrace is within the raw execution. The variation of the execution time using LiTL, in most cases, is within 2 %. Furthermore, these experiments demonstrate -especially the MG, LU, and CG kernels that generate a large number of events per second -that the instrumentation of the event-intensive applications is not harmful for their performance. The second table, Tab. IV, shows another performance metric that we used -the size of the generated traces. For the kernels that produce few events (such as FT, IS or EP), the size of traces recorded by LiTL is much smaller than of the FxT traces. The main reason for that is the header size that is relatively small for the LiTL traces. In case of the larger traces as the ones produced for LU, BT, SP or CG, the size of the header is negligible. For these kernels, the trace produced with FxT are up to 50 % larger than the LiTL traces. LiTL generates smaller traces because of the way it stores the events parameters: while FxT uses 8 bytes per parameter, LiTL only allocates the required number of bytes. The amount of the saved disk space thus depends on the number of events parameters and their type.

F. NPB-MZ

In addition, we conducted a set of experiments using the NPB-MZ kernels that combine OpenMP and MPI programming models. The combination of these programming models often occurs in practice in order to enhance the application performance and efficiently utilize multi-core nodes on standard clusters.

Tab. V presents the results of NPB-MZ-MPI benchmarks on the cluster composed of 16 Harpertown nodes. One MPI process runs on each node and 8 OpenMP threads exploit the available cores within a node. As for the NPB kernels, the EZTrace instrumentation in conjunction with both FxT and LiTL results in small performance overhead. LiTL delivers nearly the same or sometimes event better, e.g. for the SP-MZ kernel, performance as FxT. Besides, the traces produced by LiTL are much smaller in size.

VIII. CONCLUSIONS AND FUTURE WORK

In order to analyze HPC applications, there is a need for performance analysis tools with low overhead. The efficiency of such tools strongly depends on two factors: the instrumentation of functions applied and the event recording library used.

In this article, we focus on the latter issue and develop our own library named LiTL -a lightweight library for recording events -taking into account the advantages and limitations of FxT. In LiTL, we applied several optimization strategies: reducing the storage of events and traces; improving the scalability when a large number of threads record events concurrently; enabling performance analysis of parallel applications that produce hundreds or thousands of traces. These mechanisms are implemented in the LiTL library that is integrated in the EZTrace framework for performance analysis.

To evaluate the efficiency of LiTL, we conduct a set of experiments using synthetic and standard benchmarks. We compare the LiTL performance results with the ones from FxT; both event recording libraries are used within the EZTrace framework. For all test cases, LiTL either perform as good as FxT leading to low overhead or often show better results. For instance, with LiTL the size of traces is reduced up to 68 % and it also saves up to 88 % of the CPU time, especially, on multi-threaded applications. Our experiments were conducted using LiTL within EZTrace, however it is a general event recording library that can be used in conjunction with any other performance analysis tool.

Our future focus will be on integrating the litl_merge utility into the event recording process in order to produce automatically an archive of traces during the execution of the application. To reduce even more the size of the generated traces, we plan to add a compression mechanism to LiTL. Also, we aim at designing EZTrace plugins for programming models for accelerators such as CUDA and OpenCL. As new energy efficient architectures appeared, e.g. Intel Phi and Kalray Multi-Purpose Processor Array (MPPA), we plan to port both EZTrace and LiTL there.

Fig. 1 :

 1 Fig.1: Storage of different kinds of events in the trace file recorded by EZTrace. In the figure, time is the time when the event occurred; cd means the event code; tp is the event type; np stands for the number of event's parameters; par -an array of parameters.

Fig. 2 :

 2 Fig. 2: Event recording mechanism on multi-threaded applications.

Fig. 4 :

 4 Fig.4: The structure of an archive composed of multiple trace files. In the figure, fid stands for the trace file name; size is the size of a merged trace file.

TABLE I :

 I Information regarding trace files recorded by EZTrace with both FxT and LiTL on the Harpertown nodes. In the table, type is the data type of event's parameters; nb_param stands for the number of parameters.

	Type	nb_param	FxT (MB)	LiTL (MB)	Overhead (%)
	char	0	22.90	16.21	41.22
	char	1	30.53	17.17	77.81
	char	2	38.16	18.12	110.57
	char	3	45.78	19.07	140.03
	char	4	53.41	20.03	166.69
	char	5	61.04	20.98	190.94
	char	6	68.67	21.94	213.07
	int	0	22.90	16.21	41.22
	int	1	30.53	20.03	52.41
	int	2	38.16	23.84	60.03
	int	3	45.78	27.66	65.54
	int	4	53.41	31.47	69.72
	int	5	61.04	35.29	72.99
	int	6	68.67	39.10	75.63
	double	0	22.90	16.21	41.22
	double	1	30.53	23.84	28.03
	double	2	38.16	31.47	21.24
	double	3	45.78	39.10	17.09
	double	4	53.41	46.73	14.30
	double	5	61.04	54.36	12.29
	double	6	68.67	61.99	10.78

TABLE II :

 II The average time for recording one event on the Harpertown nodes. In the table, nb_th and nb_event stand for the number of threads and events, accordingly; time is in nanoseconds; VT is a shortcut for VampirTrace.

	nb_th	FxT	VT	LiTL	nb_event	Overhead (%): FxT	VT
	1	217.40	391.19	218.19	200010	-0.36	79.29
	2	488.06	395.25	220.51	400014	121.33	79.24
	3	684.71	403.28	304.58	600018	124.80	32.41
	4	946.51	417.12	380.10	800022	149.01	9.74
	5	1385.14	418.58	397.81	1000026	248.20	5.22
	6	1740.97	435.25	417.57	1200030	316.93	4.24
	7	2367.40	448.07	415.44	1400034	469.85	7.85
	8	3470.53	440.31	411.05	1600038	744.32	7.12

TABLE III :

 III The time measurements of NPB-MPI running for CLASS C problems with 64 computing processes on the Harpertown nodes. In the table, the time is in seconds.

	Bench	NPB	FxT	LiTL	Events / sec	Overhead (%): FxT	LiTL
	BT	53.75	51.29	53.31	65689	-4.58	-0.82
	MG	7.23	7.06	7.28	185491	-2.35	0.69
	LU	46.07	44.73	47.03	810584	-2.91	2.08
	FT	22.78	22.66	22.92	3681	-0.53	0.61
	IS	2.04	1.97	2.01	36619	-3.43	-1.47
	SP	80.31	77.13	76.98	58756	-3.96	-4.15
	CG	29.47	27.29	29.40	223701	-7.40	-0.24
	EP	4.99	4.98	5.07	1286	-0.20	1.60

TABLE IV :

 IV The results tracing NPB-MPI running for CLASS C problems with 64 computing processes on the Harpertown nodes.

	Bench	Events / sec	FxT (MB)	LiTL (MB)	Size difference (%)
	BT	68952	152.66	112.86	35.27
	MG	189542	45.68	32.15	42.08
	LU	838890	1866.43	1244.55	49.97
	FT	3820	0.98	0.36	171.32
	IS	35384	0.80	0.24	228.56
	SP	60699	369.22	267.87	37.84
	CG	236682	283.63	201.40	40.83
	EP	1405	0.61	0.11	439.40

TABLE V :

 V The time measurements of NPB-MZ-MPI running for CLASS C problems with 16 MPI processes and 8 OpenMP threads on the Harpertown nodes. In the table, the time is in seconds.

	Bench	NPB	FxT	LiTL	Events / sec	Overhead (%): FxT	LiTL
	LU-MZ	16.34	16.40	16.14	104836	0.37	-1.22
	BT-MZ	19.88	19.50	20.10	1252744	-1.91	1.11
	SP-MZ	22.47	22.96	22.85	2807677	2.18	1.69

This scenario represents the situation when each event's parameter can be stored using the smallest data type, which is uint8_t on X86_64 architectures.

ACKNOWLEDGMENTS

We would like to thank the AICES Graduate School and the Center for Computing and Communication of RWTH Aachen University for the provided resources.


Mehmet Ersoy
email: mehmet.ersoy@univ-tln.fr

Frederic Golay
email: frederic.golay@univ-tln.fr

Lyudmyla Yushchenko
email: lyudmyla.yushchenko@univ-tln.fr

Frédéric Golay

Adaptive scheme based on

Keywords: hyperbolic system, finite volume scheme, local mesh refinement, numerical density of entropy production, local time stepping, low density, vacuum state, strong shock, strong contact discontinuities, strong rarefaction

Introduction

We are interested in the computation of numerical entropy solutions containing vacuum states or low densities states for hyperbolic conservation laws:

∂w ∂t + ∂f (w) ∂x = 0, (t, x) ∈ R + × R w(0, x) = w 0 (x), x ∈ R. (1)
where w : R + × R → R d stands for the vector state and f : R d → R d the flux function.

It is well-known that solutions of System (1) can and will breakdown at a finite time, even if the initial data are smooth, and develop complex structure (shock wave interactions). In such a situation, the uniqueness of the (weak) solution is lost and is recovered by completing System [START_REF] Berger | Adaptive mesh refinement for hyperbolic partial differential equations[END_REF] with an entropy inequality of the form:

∂s(w) ∂t + ∂ψ(w) ∂x ≤ 0 , (2)
where (s, ψ) stands for a convex entropy-entropy flux pair. This inequality allows to select the physical relevant solution. Moreover, the entropy satisfies a conservation equation only in regions where the solution is smooth and an inequality when the solution develops shocks. In simple cases, it can be proved that the term missing in [START_REF] Chen | Vacuum states and global stability of rarefaction waves for compressible flow[END_REF] to make it an equality is a Dirac mass.

In particular, solving Equation (1) in presence of low density state or vacuum state or in presence of a strong shock also constitute a source of difficulties in regard to uniqueness of solutions to the Riemann problem. The vacuum states introduces an indeterminacy in the variables which describe the flow (vacuum has zero energy, momentum and pressure and the velocity is not defined since u = q/ρ where q stands for the momentum) while low density state often yields to numerical instabilities. Strong shocks induce a large discontinuity in the pressure and often one has to impose limitations on the strength of shocks from mathematical as well as numerical viewpoint (see e.g. [START_REF] Landau | Fluid mechanics[END_REF][START_REF] Sachdev | Chapman & Hall/CRC Monographs and Surveys in Pure and Applied Mathematics[END_REF][START_REF] Dafermos | Hyperbolic conservation laws in continuum physics[END_REF]). In general, it often yields singularities in the physical systems producing mathematical ill-posedness and numerical instabilities for solutions which are close to those states. We refer to [START_REF] Serre | Systems of conservation laws[END_REF][START_REF] Chen | Vacuum states and global stability of rarefaction waves for compressible flow[END_REF][START_REF] Huang | Convergence rate for compressible Euler equations with damping and vacuum[END_REF][START_REF]Asymptotic behavior of the solutions to the damped compressible Euler equations with vacuum[END_REF][START_REF] Leveque | Finite volume methods for hyperbolic problems[END_REF][START_REF] Dafermos | Hyperbolic conservation laws in continuum physics[END_REF] for mathematical discussions about those states. As a consequence, to compute accurately such type of solutions is a challenging problem that we intend to implement efficiently in a multi-scale adaptive framework.

Following Ersoy et al [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF] and the above considerations, we will numerically investigate the robustness of the adaptive numerical scheme through the computation of those solutions. Let us recall that the adaptive multi-scale scheme is based on the numerical density of the entropy production. This quantity is a measure of the amount of violation of the entropy equation (as a measure of the local residual as already pointed out in [START_REF] Berger | Adaptive mesh refinement for hyperbolic partial differential equations[END_REF][START_REF] Houston | A posteriori error analysis for numerical approximations of Friedrichs systems[END_REF][START_REF] Karni | A smoothness indicator for adaptive algorithms for hyperbolic systems[END_REF][START_REF] Karni | Local error analysis for approximate solutions of hyperbolic conservation laws[END_REF][START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF]). Thus, the numerical density of entropy production provides information on the need to locally refine the mesh (e.g. if the solution develops discontinuities) or to coarsen the mesh (e.g. if the solution is smooth and well-approximated). We have shown that for solutions far from those states, such indicator is able (see Ersoy et al [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF]) not only to provide an efficient a posteriori error, but also to reproduce the qualitative structure of the solution and to pilot the adaptive scheme. Let us also emphasize that we prefer to use the classical Godunov solver instead of linearized or approximated solver in the finite volume framework. Indeed, it may happen for data close to the vacuum state or in situations where there is a strong expansion, those linearized scheme fail (see for instance, citeEMRS91,leveque2002,Toro09).

The paper is organized as follows. In Section 2, we fix notations and we recall the one dimensional gas dynamics equations for perfect gas. In section (3), we recall the finite volume approximation in the multi-scale adaptive framework as introduced by Ersoy et al [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF]. Moreover, we recall the procedure to deal with solutions close to vacuum states. Finally, in Section 4, through several numerical test case, we show the efficiency of the present numerical scheme for solutions which are close to vacuum state or/and in presence of strong discontinuities and rarefactions. In particular, we numerically show that those solutions are very well-computed. To illustrate this statement, we study the qualitative behavior of the adaptive numerical scheme with respect to the L 1 error norm providing quantitative results. We only restrict ourselves to solutions which are only C 0 and even with strong discontinuities as for shock or contact discontinuities since for smooth solutions we recover the order p for a numerical scheme of order p. Even in presence of severe test cases, we emphasize that the present scheme is super-convergent or at least of order 1. We focus on the following tests:

• a modified Sod test problem,

• a 123 problem,

• a left half of the blast wave of Woodward and Colella,

• a Two shock Riemann problem,

• a modified Leblanc test case,

• a blast wave problem of Woodward and Colella to assess the ability to conserve the entropy inequality at the discrete level, the performance of the scheme for low density flows and the robustness of the numerical method. One can found other test cases for the one dimensional gas dynamics equation for perfect gas, for instance, in [START_REF] Woodward | Trade-offs in designing explicit hydrodynamical schemes for vector computers[END_REF][START_REF]The numerical simulation of two-dimensional fluid flow with strong shocks[END_REF][START_REF] Einfeldt | On Godunov-type methods near low densities[END_REF][START_REF] Gallouët | Some recent finite volume schemes to compute Euler equations using real gas EOS[END_REF][START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics[END_REF].

The 1D gas dynamics equations

Let us consider a compressible perfect fluid confined in a domain ⊗. The governing equations for the motion of the compressible fluid in [0, T ] × ⊗ are the so called Euler equations

             ∂ρ ∂t + ∂ρu ∂x = 0, ∂ρu ∂t + ∂ ρu 2 + p ∂x = 0, ∂ρE ∂t + ∂ (ρE + p) u ∂x = 0 (3)
where ρ is the density, p = (γ -1)ρε is the pressure, u the velocity and E the total energy defined by:

E = ε + u 2 2 .
Here, ε stands for the internal specific energy and γ is the ratio of the specific heats. Making use of the notations

w =   ρ ρu ρE   , f (w) =   ρu ρu 2 + p (ρE + p) u  
for the conservative variables and the conservative flux, the system (3) is written as a system of conservation laws (1). We complete System (3) with the entropy inequality [START_REF] Chen | Vacuum states and global stability of rarefaction waves for compressible flow[END_REF] where the convex entropy s(w) and the entropy flux pair ψ(w) are classically given by the following relations:

s(w) = -ρ ln p ρ γ , ψ(w) = u s(w) . (4)
Let us note that even if System (3) is strictly hyperbolic on the set {ρ > 0}, the previous quantities (4) make a sense for ρ 0.

The adaptive numerical scheme

For the sake of completeness, we recall in this section, the first order adaptive multi scale scheme based on numerical density of entropy production for conservation laws introduced by the authors in [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF]. Interested reader can found the details of the second order scheme with MUSCL reconstruction as well as the second order Adams-Bashforth time integration in [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF].

One-dimensional finite volume formulation.

We recall here the well-known construction of a numerical approximation of the general non linear hyperbolic conservation laws [START_REF] Berger | Adaptive mesh refinement for hyperbolic partial differential equations[END_REF] The computational domain is split into control volumes

C k = x k-1/2 , x k+1/2 of mesh size h k with x k±1/2 = x k ± h k /2.
The unknowns w(t, x) are approximated by their mean values on the cell C k at time t:

w k (t) ≃ 1 h k C k w (t, x) dx .
Noting, t n+1 = t n + δt n the discrete time and Integrating (3) over each strip (t n , t n+1 ) × C k and applying the Green's formula we obtain:

w n+1 k -w n k δt n + δF n k h k = 0 , (5)
with

δF n k := F n k+1/2 -F n k-1/2
where F n k±1/2 stands for the numerical fluxes at the interface x k±1/2 . We use here the Godunov solver, i.e. we defined

F n k±1/2 (w n k , w n k+1 ) = f R(0 ± , w n k , w n k+1 )
where f (R(0 ± , w L , w R )) is determined from the exact solution R(0 ± , w L , w R ) of the Riemann problem associated with the left w L and the right w R state:

     ∂w ∂t + ∂f (w) ∂x = 0 , w 0 (x) = w L if x < 0 , w R if x > 0 . (6)
In the same way, the discretization of the entropy inequality (4) yields to:

S n k := s(w n+1 k ) -s(w n k ) δt n + δψ n k h k (7)
where

δψ n k := ψ n k+1/2 -ψ n k-1/2 := ψ R(0 -, w n k , w n k+1 ) -ψ R(0 + , w n k-1 , w n k ) . The quantity S n
k defines the so-called numerical density of entropy production and the value

P = n,k S n k δt n h k . ( 8
)
is the total numerical entropy production.

Mesh refinement.

In order to compute accurately the solution of hyperbolic systems, the strategy [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF] consists to use the numerical density of entropy production as a mesh refinement parameter in the following meaning.

We present here a local mesh refinement procedure driven by the numerical density of entropy production. In order to reduce the time necessary to manage the refinement, we use "macro cells" which could be refined by generating hierarchical grids. Each cell can be split in two sub-cells. We thus produce a dyadic cells graph, whose numbering (in basis 2) allows a quick computing scan to determine the adjacent cells. For the convenience of the reader, we make use of the following notations: let k b be the index which makes reference to the macro cell numbered k and b a binary number which contains the hierarchical information of a sub-cell. In particular, the level of a sub-cell C k b is defined as the length(b) -1. For instance, a macro cell C k0 of level 0 will be split into two sub-cells C k00 and C k01 of level 1. A mesh refinement example is proposed in Figure 1. The mesh refinement procedure can be simply expressed as follows. We first set a mesh refinement parameter S. For instance, it can be the mean value over the domain Ω at time t n :

S = 1 |Ω| k b S n k b . (9)
We then define two coefficients 0 α min α max 1, which determine the ratio of numerical production of entropy leading to mesh refinement or mesh coarsening. Thus, for each cell C k b :

• if S n k b > Sα max , the mesh is refined and split into two sub-cells C k b0 and C k b1 ,

• if S n k b0 < Sα min and S n k b1 < Sα min , the mesh is coarsened into a cell C k b .
On one hand, if a cell C k b is split into two sub-cells C k b0 and C k b1 , averaged values at time t n are projected on each sub-cell (for the first order scheme):

w n k b0 = w n k b1 = w n k b
and the numerical fluxes between the two sub-cells are defined as:

           F n k b0 -1/2 = F n k b -1/2 , F n k b0 +1/2 = F n k b1 -1/2 = f (w n k b ) , F n k b1 +1/2 = F n k b +1/2 .
The adaptive mesh refinement scheme is also written in a local time stepping method to gain in time (we refer to [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF] for more details).

Numerical results

In this section, we show the efficiency of the adaptive numerical1 scheme based on the numerical density of entropy production through severe numerical experiments. In [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF], we have pointed out that the use of the numerical density entropy production as a mesh refinement parameter is a relevant local error indicator (everywhere where the solution remains smooth) and a discontinuity detector (shocks and oscillating solutions are very well-captured). Far from large discontinuities and vacuum states, we have performed several numerical tests (the Sod's shock tube problem, the Lax's shock tube problem and the Shu-Osher test problem) and we have shown an impressive improvement with respect to uniform grids even if a large number of uniform cells are used. Moreover, the adaptive strategy, combined with a local time stepping algorithm, can significantly reduce the computational time keeping the same order of accuracy. We also refer to the 2D and 3D extension of the work [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF] in [START_REF] Yushchenko | Entropy production and mesh refinement -Application to wave breaking[END_REF].

From now on, we focus on severe numerical experiments designed to assess the performance of the numerical method in presence of strong discontinuities (both shocks and contact discontinuities), strong rarefactions as well as for low density flows (see for instance [START_REF] Sod | A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws[END_REF][START_REF] Woodward | High resolution difference schemes for compressible gas dynamics[END_REF][START_REF]The numerical simulation of two-dimensional fluid flow with strong shocks[END_REF][START_REF] Einfeldt | On Godunov-type methods near low densities[END_REF][START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics[END_REF]). Therefore, in the sequel, we focus on the robustness of the present adaptive scheme.

For all test cases, numerical solutions are computed using the one-dimensional gas dynamics equations for ideal gas (3) and we make use of the following notations and settings:

Notations and settings

• In what follows, we perform several numerical tests using first and second order schemes. Thus, we will refer to AB1 as the first order scheme, AB2 as the second order Adams-Basforth scheme, RK2 as the second order Runge-Kutta scheme. AB2 and RK2 use a MUSCL reconstruction. Moreover, all computations are made with a dynamic grid except if the acronym ends with the capital letter "U" which refers to a uniform fixed grid. All computations are made with a uniform time step except if the acronym ends with the capital letter "M" which refers to the local time stepping algorithm (see [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF] for further details on the local time stepping adaptation). We also note the cpu time cpu g for the uniform (or global) time stepping, -cpu l for the local time stepping .

• We will compare the adaptive numerical solution to the one computed on a uniform fixed grid. To have a coherent support of comparison, the solution on the fixed grid will be computed with N Lmax cells. N Lmax stands for the average number of cells used during a simulation of an adaptive scheme with a maximum level L max .

• All of the presented results display the density which is positive. Thus, for the convenience of the reader, in all figures the numerical density of entropy production is plotted with the reversed sign to overlap with the density.

• For all numerical tests, the initial mesh is adapted to the initial data during few iterations (except for the numerical solution computed on uniform fixed grid) and we have used the following threshold parameters: Mesh refinement parameter α max : 0.01 , Mesh coarsening parameter α min : 0.001 ,

Mesh refinement parameter S : 1 |Ω| k b S n k b .
• For each test 1 -test 5, X ex will make reference to the exact solution of the Riemann problem computed using 5 000 uniform cells. For test 6, we use the notation X ref which stands for the reference solution computed on a fixed grid with 20 000 cells using the RK2 scheme.

• To study the numerical convergence of AB1,AB1U, AB1M, RK2U, AB2U, AB2M and RK2 schemes, we have used the discrete l 1

x -norm of the density, momentum, pressure and the internal energy error. We generally used the convergence rate of the internal energy which is more representative of the efficiency of the scheme especially for low density flows.

• On the tables of performance 1, 4, 7, 10, 13 and 16, since the error of the AB1 and AB1M schemes are similar, we only use the results of the AB1 scheme. We also use the AB2 scheme results for the table of performance 2, 5, 8, 11, 14 and 17.

Test 1: a modified version of the Sod's shock tube problem

We first consider a modified version of the popular Sod's test [START_REF] Sod | A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws[END_REF] (see also [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF]) which models the onedimensional flow resulting from the rupture of a membrane separating air with different density/pressure inside a shock tube with a null speed. Besides the fact that the exact solution can be easily computed, the solution of the modified Sod's shock tube problem consists of a right shock wave, a right traveling contact wave and a left sonic rarefaction wave. This is useful to test the ability of the scheme to preserve in-cell entropy inequality even if the adaptive strategy is used. The test is given by a one-dimensional Riemann problem with the following data:

x ∈ [-1, 1], (ρ, u, p)(0, x) = (1, 0.75, 1), x 0, (0.125, 0, 0.1), x > 0 .

The computational domain here is [-1, 1] with prescribed free boundary conditions. For each numerical computation, the following parameters have been used:

CFL : 0.25, Simulation time (s) : 0.2,

Initial number of cells

: 200, Maximum level of mesh refinement : L max .

In figure 2, we display the solution profiles (as well as the exact one) for the density (figure 2 As usual, by the use of the Godunov solver, the shock wave is rather well-computed with a smearing effects over neighborhood cells for the solution on a fixed grid while it perfectly match when the adaptive grid technique is used (as observed for each state (ρ, u, p, ε) on figure 2). The numerical diffusion decreases when the solution is refined in around the shock wave as predicted by the numerical density of entropy production, see figure 2(a). Let us also add that no spurious oscillations are observed in the vicinity of the shock when the adaptive grid strategy is used. For more details we refer to [5, Section 5.1.1.].

The contact discontinuity, which is much more difficult to resolve accurately, is here efficiently computed with an adaptive grid since the mesh is well-refined. Another positive feature of the numerical approximation is that the speed of the contact discontinuity is correct for the solution on an adaptive grid even if we use the Osher-Sanders projections which is well-known to be locally non conservative and non consistent (see [START_REF] Osher | Numerical approximations to nonlinear conservation laws with locally varying time and space grids[END_REF] and [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF], since the correct speed of propagation of contact discontinuities are well-computed when the scheme is conservative.)

Another feature concerns the rarefaction waves which are smooth flows. In general, any numerical Finite Volume schemes are able to well-approximate this smooth solution except near the head and the tail. Since the numerical density of entropy production is non zero near the head and the tail (see figure 2(a)), the adaptive solution corresponds to a-well-approximation of the exact solution. Moreover, let us note that the so-called entropy glitch (discontinuity inside the rarefaction) is observed only on a fixed grid (see figures 2(a), 2(b) and 2(d)).

In the following, we display on figure 3 the numerical order of convergence using the discrete l 1

x -norm on the error of the internal energy for the first and second order schemes. We also complete the previous results with the l 1

x error of ρ, p, u and ε on table 1 and 2 computed with the first and the second order scheme. We finally summarize on table 3 all the rate of convergence of all numerical schemes on ρ, p, u and ε. x with respect to the average number of cells at time t = 0.2. x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time. As expected, since the solution is discontinuous, we recover approximately order 1 for the first and second order method on uniform mesh. Moreover, we also get "super-convergence" for first order adaptive schemes (see figure 3(a) and tables 1 and 2) as already obtained in [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF] for the internal energy. The computational time is considerably reduced when the local time stepping algorithm is used as observed on table 1 and table 4.

-2.2 -2 -1.8 -1.6 -1.4 -1.2 -1 2.
L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.163E-
L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.116E-

Test 2: a modified 123 problem

The 123 numerical experiment or the Two-Rarefaction test problem (see for instance [START_REF] Einfeldt | On Godunov-type methods near low densities[END_REF]) is designed to test the ability of the scheme for low density flows. It consists in two symmetric strong rarefaction waves and a trivial contact wave of zero speed. In this test case, the pressure in the star region (i.e., the between the non-linear waves) is close to zero. As a consequence, it yields to some technical difficulty to find the pressure in the iterative Newton-Raphson procedure if the initial guess is not appropriately chosen. In that case, we have used the guess resulting from a linearized solution based on primitive variables (see [START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics[END_REF]Formula 4.47]):

p 0 = max(α, P ) where P = 1 2 (p l + p r ) - 1 8 (u l -u r ) (ρ l + ρ r ) (c l + c r )
where X l (resp. X r ) stands for the left (resp. the right) state, c = (γ -1)p ρ for the sonic speed and α is a given tolerance (10 -12 in the present case). Instead of using the original test,

x ∈ [-1, 1], (ρ, u, p)(0, x) = (1, -2, 0.4), x 0, (1, 2, 0.4), x > 0

we consider the modified one given by the following Riemann data:

x ∈ [-1, 1], (ρ, u, p)(0, x) = (1, -3, 0.3), x 0, (1, 3, 0.3), x > 0
where the densities is now close to 10 -3 (instead of 10 -1 in the original case). The computational domain here is [-1, 1] with free boundary conditions and for each numerical computation, the following parameters have been used: CFL : 0.25, Simulation time (s) : 0.15, Initial number of cells : 200, Maximum level of mesh refinement : L max .

In figure 4, we display the solution profiles (as well as the exact one) for the density (figure 4 As pointed out in [START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics[END_REF], the Godunov's solver is generally accurate as regards of the primitive variables and even better with the adaptive strategy. Nevertheless, it is no more true when we look for the specific internal energy, ε = p (γ -1)ρ , as one can observe on figure 4(d). We also plot the solution given by the RK2 scheme (which is well-known to be accurate) using 5 000 uniform cells on figure 4(d).

Even in presence of strong rarefactions (see figures 4(a) and 4(b)), the numerical density of entropy production is non zero and thus as pointed out in the previous section, except the head and the tail of the rarefaction, the solution is rather well-computed and even more in the adaptive framework (both from ρ, ρu, p and ε).

Let us also emphasize that in the present adaptive scheme, besides the particular guess in the iterative scheme for finding p, we do not make specific treatment when the solution is close to the vacuum. In particular, in the low density area, the numerical density of entropy production almost vanishes as shown in figure 4(a).

As done in the previous section, we now complete the previous results with table of performances 4-5 and 6 as well as the rate of convergence of the first and the second order schemes for the internal energy, see figure 5. x with respect to the average number of cells at time t = 0.15. x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time. As expected, the numerical order of the adaptive scheme are better but we do not observe a specific "super-convergence" since no discontinuities are available in the density. Let us also add that the behavior of the internal energy (see figure 4(d)) can be explained by the fact that the density converges faster than the pressure to 0.

-1.4 -1.3 -1.2 -1.1 -1 -0.
L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.451E-
L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.357E-

Test 3: the left half of the blast wave of Woodward and Colella

We now focus on a very severe test problem. This test is designed to assess the efficiency of the numerical scheme when the solution contains a strong right shock, a left rarefaction and a contact surface. One can also regard this test problem as the left half of the blast wave problem of Woodward and Colella [START_REF]The numerical simulation of two-dimensional fluid flow with strong shocks[END_REF] (see also Test 6).

The test 3 is given by the use of the following Riemann data:

x ∈ [-1, 1], (ρ, u, p)(0, x) = (1, 0, 1000), x 0, (1, 0, 0.01), x > 0 .

The computational domain here is [-1, 1] and we use free boundary conditions. For each numerical computation, we use the following parameters: CFL : 0.25, Simulation time (s) : 0.02, Initial number of cells : 200, Maximum level of mesh refinement : L max .

In figure 6, we display the solution profiles (as well as the exact one) for the density (figure 6(a)), the pressure (figure 6(b)), the momentum (figure 6(c)) and the internal energy (figure 6(d)) of the numerical solution computed with the first order scheme AB1U and AB1 computed with L max = 5. The numerical density of entropy production (S n k ) is also plotted on figure 6(a). Even if the solution contains a non physical large shock, this severe test case allows to numerically show the robustness of the numerical scheme. In particular, it is difficult to capture efficiently the pressure (as well as the density) due to the large amplitude of the shock. To get an accurate approximation of the exact solution, a reasonable number of cells (N Lmax = 353) is required to well-approximate the post shock as clearly seen on the density figure 6(a), the pressure figure 6(b) and the internal energy figure 6(d) (and even better for second order scheme, see table 7 and8). As observed on figure 6(d), we expect a "super-convergence" (as in the first test problem). We refer to tables 7, 8 and 9. x with respect to the average number of cells at time t = 0.02. x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.329E+00 0.
L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.252E+00

Test 4: the Two Shock Riemann problem

As the test 3, test 4 is also a severe test. The solution is composed of three strong discontinuities (two shocks and a contact discontinuity) traveling to the right. The left shock has a small positive propagation speed while the right is the fastest wave. The Riemann data are here:

x ∈ [-1, 1], (ρ, u, p)(0, x) = (5.99924, 19.5975, 460.894), x 0, (5.99242, -6.19633, 46.095), x > 0 .

The computational domain here is [-1, 1] and we use free boundary conditions. For each numerical computation, we use the following parameters: CFL : 0.25, Simulation time (s) : 0.035, Initial number of cells : 200, Maximum level of mesh refinement : L max .

In figure 8, we display the solution profiles (as well as the exact one) for the density (figure 8(a)), the pressure (figure 8(b)), the momentum (figure 8(c)) and the internal energy (figure 8(d)) of the numerical solution computed with the first order scheme AB1U and AB1 computed with L max = 5. The numerical density of entropy production (S n k ) is also plotted on figure 8(a). Even if the left shock wave is very well-resolved, the right one (which is the fastest wave) is smeared over several cells while the numerical diffusion for the contact discontinuity is important. The smearing effects can be reduced with the adaptive scheme since in that case the numerical density entropy production detects the presence of huge discontinuities (see figures 8(a) and 8(d)) and therefore the mesh is refined. Let us also note the presence of spurious oscillations in the slowly moving left shock which are considerably reduced in the adaptive framework (see figure 8(c)). As in the previous test case, the rate of convergence (see table 12) as well as the table of performance (see table 10 x with respect to the average number of cells at time t = 0.035. x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.103E+01 0.
L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.732E+00

Test 5: a modified Leblanc test problem

We present in this section a modified version of the so-called Leblanc test problem which is an extreme shock tube problem. The initial discontinuity separates a region of very high energy and density from one of low energy and density. This is a much more severe test than the Sod Problem. The solution consists of a rarefaction moving to the left, and a contact discontinuity and a strong shock moving to the right. The Riemann data for this test problem are:

x ∈ [-1, 1], (ρ, u, p)(0, x) = (1, 0, 0.1), x 0, (0.001, 0, 10 -9 ), x > 0 .

The computational domain here is [-1, 1] with prescribed free boundary conditions and we have used the following parameters: CFL : 0.25, Simulation time (s) : 0.7, Initial number of cells : 200, Maximum level of mesh refinement : L max .

In figure 10, we display the solution profiles (as well as the exact one) for the density (figure 10 This modified test is a challenging test problem since the presence of the strong left rarefaction wave and discontinuities of small amplitude make it difficult. Indeed, while the speed of the contact discontinuity is correct, we make an error on the correct shock speed as one can observe on figures 10(c) and 10(d) when the number of cells is approximately 1 000 cells. As a consequence, the post shock values are extrapolated (see, 10(d)) even if the contact discontinuity and the speed of propagation of the rarefaction are well-computed.

Starting with 500 cells instead of 200 cells, we get better results as displayed on figure 11(d). This behavior is not due to a lack of accuracy of the present scheme but it is due to the extreme difficulty of this test problem. Let us note that in the literature for the original Leblanc test problem, they generally use more than 3 000 cells with high order scheme to resolve it "accurately". x with respect to the average number of cells at time t = 0.15.

L max ρ -ρ ex l 1 x p -p ex l 1 x u -u ex l 1 x ε -ε ex l 1 x N T f cpu g cpu l 1 0.151E-

Test 6: the blast wave problem

The Woodward-Colella blast wave test problem can be extremely difficult to solve on a uniform grid even if we use a very large number of cells. This one-dimensional test problem was initially introduced in [START_REF] Woodward | Trade-offs in designing explicit hydrodynamical schemes for vector computers[END_REF] and it makes one of the main difficult and the main challenging test case to solve in order to show the efficiency of the numerical method. In particular, it illustrate the strong relationship between the accuracy of the overall flow solution and the thinness of discontinuities on the grid involving multiple interactions of discontinuities (shocks and contact discontinuities) and rarefactions with each other. The initial condition consists of three constant states

x ∈ [-1, 1], ρ(0, x) = 1, u(0, x) = 0, p(0, x) =    1000, x 0.1, 0.01, 0.1 < x 0.9, 100, x > 0.9 .

on a computational domain [0, 1] with prescribed reflecting boundary conditions. In that case, as a reference solution, we compute the solution on a uniform grid with 20 000 cells using the RK2U scheme. We also use the following parameters: CFL : 0.25, Simulation time (s) : 0.038, Initial number of cells : 200, Maximum level of mesh refinement : L max .

In figure 13, we display the solution profiles (as well as the reference one) for the density (figure 13(a)), the pressure (figure 13(b)), the momentum (figure 13(c)) and the internal energy (figure 13(d)) of the numerical solution computed with the first order scheme AB1U and AB1 computed with L max = 5. The numerical density of entropy production (S n k ) is also plotted on figure 13(a). x with respect to the average number of cells at time t = 0.038.

L max ρ -ρ ex l We numerically show that the adaptive grid strategy using the numerical density of entropy production can considerably increase the accuracy without requiring a large number of cells and without a lot of extra CPU times as required by solving the problem on a uniform grid.

Concluding remarks and perspectives

In this paper, we have used the adaptive numerical scheme based on the numerical density of entropy production introduced by Ersoy et al [START_REF] Ersoy | Adaptive multiscale scheme based on numerical density of entropy production for conservation laws[END_REF]. This scheme employs a useful adaptive strategy based on the mathematical (physical) entropy which indicates on the need to refine or to coarsen the mesh, in particular in the presence of strong discontinuities. Moreover, combining into a local time stepping method, one can considerably decrease the computational cost keeping the same order of accuracy. Through several severe test problems for which the solution contains strong discontinuities, strong rarefaction and/or low density, we have numerically shown the robustness of the scheme. Indeed, resulting numerical solutions are both accurate for discontinuities and the head/tail of rarefaction waves which are well-known to be difficult to handle numerically.

Figure 1 :

 1 Figure 1: Example of hierarchical dyadic tree.

  (a)), the pressure (figure2(b)), the momentum (figure2(c)) and the internal energy (figure2(d)) of the numerical solution computed with the first order scheme AB1U and AB1 computed with L max = 5. The numerical density of entropy production (S n k ) is also plotted on figure2(a).ρu on adaptive mesh with L max = 5 (N L max = 449) ρu on fixed mesh N = 449 (ρu) ex (c) Momentum.

ε

  on adaptive mesh with L max = 5 (N L max = 449) ε on fixed mesh N = 449 ε ex (d) Internal energy.

Figure 2 :

 2 Figure 2: Test 1: numerical solution at time t = 0.2 s with L max = 5 computed with the AB1 scheme.

Figure 3 :

 3 Figure 3: Test 1: ε -ε ex l 1x with respect to the average number of cells at time t = 0.2.

ε

  (a)), the pressure (figure 4(b)), the momentum (figure 4(c)) and the internal energy (figure 4(d)) of the numerical solution computed with the first order scheme AB1U and AB1 computed with L max = 5. The numerical density of entropy production (S n k ) is also plotted on figure 4(a). ρu on adaptive mesh with L max = 5 (N L max = 1084) ρu on fixed mesh N = 1084 ρu ex on adaptive mesh with L max = 5 (N L max = 1084) ε on fixed mesh N = 1084 ε on fixed mesh with RK2 N = 5000 ε ex (d) Internal energy.

Figure 4 :

 4 Figure 4: Test 2: numerical solution at time t = 0.15 s with L max = 5 computed with the AB1 scheme.

Figure 5 :

 5 Figure 5: Test 2: ε -ε ex l 1x with respect to the average number of cells at time t = 0.15.

ε

  ρu on adaptive mesh with L max = 5 (N L max = 353) ρu on fixed mesh N = 353 ρu ex on adaptive mesh with L max = 5 (N L max = 353) ε on fixed mesh N = 353 ε ex (d) Internal energy.

Figure 6 :

 6 Figure 6: Test 3: numerical solution at time t = 0.02 s with L max = 5 computed with the AB1 scheme.

Figure 7 :

 7 Figure 7: Test 3: ε -ε ex l 1x with respect to the average number of cells at time t = 0.02.

ε

  on adaptive mesh with L max = 5 (N L max = 292) ε on fixed mesh N = 292 ε ex (d) Internal energy.

Figure 8 :

 8 Figure 8: Test 4: numerical solution at time t = 0.035 s with L max = 5 computed with the AB1 scheme.

Figure 9 :

 9 Figure 9: Test 4: ε -ε ex l 1x with respect to the average number of cells at time t = 0.035.

ε

  (a)), the pressure (figure 10(b)), the momentum (figure 10(c)) and the internal energy (figure 10(d)) of the numerical solution computed with the first order scheme AB1U and AB1 computed with L max = 5. The numerical density of entropy production (S n k ) is also plotted on figure 10(a). ρu on adaptive mesh with L max = 5 (N L max = 995) ρu on fixed mesh N = 995 ρu ex on adaptive mesh with L max = 5 (N L max = 995) ε on fixed mesh N = 995 ε ex (d) Internal energy.

Figure 10 :

 10 Figure 10: Test 5: numerical solution at time t = 0.7 s with L max = 5 computed with the AB1 scheme.

ρε

  on adaptative mesh with L max = 5 (N L max = 3575) ρ on fixed mesh N = 3575 ρ ex Τ (a) Density and numerical density of entropy production. ρu on adaptative mesh with L max = 5 (N L max = 3575) ρu on fixed mesh N = 3575 ρu ex on adaptative mesh with L max = 5 (N L max = 3575) ε on fixed mesh N = 3575 ε ex (d) Internal energy.

Figure 11 :

 11 Figure 11: Test 5: numerical solution at time t = 0.7 s with L max = 5 computed with the AB1 scheme with 500 cells.

Figure 12 :

 12 Figure 12: Test 5: ε -ε ex l 1x with respect to the average number of cells at time t = 0.15.

Figure 14 :

 14 Figure 14: Test 6: ε -ε ex l 1x with respect to the average number of cells at time t = 0.038.

Table 1 :

 1

		01	0.112E-01	0.302E-01	0.813E-01	200 0.28 0.28
	2	0.103E-01	0.617E-02	0.154E-01	0.515E-01	246 0.66 0.37
	3	0.631E-02	0.331E-02	0.761E-02	0.315E-01	307 1.66 0.66
	4	0.354E-02	0.167E-02	0.345E-02	0.175E-01	399 3.18 1.36
	5	0.142E-02	0.653E-03	0.121E-02	0.679E-02	567 7.37 3.85

Test 1 : Convergence tests for first order AB1 scheme, l 1

Table 2 :

 2 Test 1 : Convergence tests for second order AB2 scheme, l 1 x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

		01	0.725E-02	0.197E-01	0.591E-01	200 0.26 0.26
	2	0.736E-02	0.396E-02	0.101E-01	0.378E-01	238 0.58 0.35
	3	0.460E-02	0.219E-02	0.509E-02	0.235E-01	291 1.26 0.64
	4	0.265E-02	0.115E-02	0.239E-02	0.134E-01	384 2.82 1.38
	5	0.116E-02	0.545E-03	0.998E-03	0.541E-02	545 7.34 4.12

Table 3 :

 3 Test 1 : Convergence tests for first and second order.

Table 4 :

 4 Test 2 : Convergence tests for first order AB1 scheme, l 1

		01	0.214E-01	0.233E+00 0.368E+00 200	0.20	0.20
	2	0.290E-01	0.137E-01	0.161E+00 0.270E+00 318	0.57	0.43
	3	0.181E-01	0.822E-02	0.125E+00 0.184E+00 546	1.67	1.26
	4	0.104E-01	0.435E-02	0.827E-01	0.114E+00 986	6.00	4.97
	5	0.403E-02	0.156E-02	0.359E-01	0.543E-01 1856 22.57 19.68

Table 5 :

 5 Test 2 : Convergence tests for second order AB2 scheme, l 1 x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

		01	0.169E-01	0.189E+00 0.311E+00 200	0.20	0.20
	2	0.231E-01	0.106E-01	0.154E+00 0.221E+00 314	0.56	0.44
	3	0.147E-01	0.630E-02	0.112E+00 0.148E+00 536	1.81	1.44
	4	0.853E-02	0.337E-02	0.701E-01	0.914E-01	980	6.40	5.38
	5	0.340E-02	0.125E-02	0.309E-01	0.449E-01 1842 22.98 19.85

Table 6 :

 6 Test 2 : Convergence tests for first and second order.

Table 7 :

 7 Test 3 : Convergence tests for first order AB1 scheme, l 1

			268E+02 0.119E+01 0.134E+03 200 0.37 0.36
	2	0.225E+00 0.201E+02 0.830E+00 0.986E+02 222 0.78 0.47
	3	0.142E+00 0.153E+02 0.590E+00 0.703E+02 253 1.74 0.75
	4	0.819E-01	0.123E+02 0.445E+00 0.461E+02 286 4.06 1.50
	5	0.345E-01	0.967E+01 0.338E+00 0.236E+02 322 9.64 3.50

Table 8 :

 8 Test 3 : Convergence tests for second order AB2 scheme, l 1 x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

			0.182E+02 0.813E+00 0.103E+03 200 0.38 0.39
	2	0.164E+00 0.131E+02 0.540E+00 0.747E+02 218 0.82 0.48
	3	0.103E+00 0.103E+02 0.390E+00 0.527E+02 240 1.75 0.74
	4	0.600E-01	0.865E+01 0.311E+00 0.344E+02 265 3.83 1.38
	5	0.261E-01	0.794E+01 0.276E+00 0.186E+02 298 8.80 2.72

Table 9 :

 9 Test 3 : Convergence tests for first and second order.

  and table11) are exceptional since only 300 cells are used.

				L 1 discretisation error					L 1 discretisation error
		1			AB1U : order = 0.57887		1			RK2U : order = 0.6143
					AB1 : order = 5.4405					AB2U : order = 0.6112
		0.8			AB1M : order = 5.0193		0.8			RK2 : order = 6.7297 AB2 : order = 6.5196
											AB2M : order = 3.9711
								0.6
		0.6
	log10(error)	0.4					log10(error)	0.2 0.4
		0.2
								0
		0						-0.2
		-0.2						-0.4
		2.3	2.35	2.4	2.45	2.5	2.55	2.3	2.35	2.4	2.45	2.5	2.55	2.6
				log10(mesh-points)						log10(mesh-points)
			(a) First order scheme.					(b) Second order scheme.

Table 10 :

 10 Test 4 : Convergence tests for first order AB1 scheme, l 1

			266E+02 0.279E+00 0.860E+01 200 0.34 0.34
	2	0.641E+00 0.134E+02 0.142E+00 0.568E+01 225 0.72 0.43
	3	0.389E+00 0.674E+01	0.687E-01	0.356E+01 249 1.61 0.66
	4	0.211E+00 0.270E+01	0.300E-01	0.199E+01 276 3.60 1.23
	5	0.794E-01	0.779E+00	0.953E-02	0.745E+00 312 7.78 3.15

Table 11 :

 11 Test 4 : Convergence tests for second order AB2 scheme, l 1 x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

			0.175E+02 0.185E+00 0.648E+01 200 0.36 0.36
	2	0.461E+00 0.929E+01	0.985E-01	0.423E+01 218 0.78 0.46
	3	0.281E+00 0.471E+01	0.470E-01	0.262E+01 237 1.72 0.69
	4	0.153E+00 0.197E+01	0.223E-01	0.148E+01 259 3.62 1.25
	5	0.579E-01	0.679E+00	0.881E-02	0.565E+00 290 7.97 4.18

Table 12 :

 12 Test 4 : Convergence tests for first and second order.

Table 13 :

 13 Test 5 : Convergence tests for first order AB1 scheme, l 1 x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

		01	0.166E-02	0.390E-01	0.413E-01	200	0.30	0.31
	2	0.964E-02	0.102E-02	0.271E-01	0.298E-01	313	0.92	0.74
	3	0.601E-02	0.599E-03	0.209E-01	0.219E-01	524	2.98	2.30
	4	0.350E-02	0.332E-03	0.157E-01	0.150E-01	912	9.68	7.87
	5	0.151E-02	0.138E-03	0.112E-01	0.878E-02 1642 32.13 26.42
	L max	ρ -ρ ex l 1 x	p -p ex l 1 x	u -u ex l 1 x	ε -ε ex l 1 x	N T f	cpu g	cpu l
	1	0.103E-01	0.106E-02	0.304E-01	0.333E-01	200	0.35	0.35
	2	0.680E-02	0.669E-03	0.240E-01	0.253E-01	311	1.03	0.81
	3	0.443E-02	0.415E-03	0.189E-01	0.189E-01	510	3.09	2.39
	4	0.277E-02	0.250E-03	0.146E-01	0.133E-01	887	9.58	7.67
	5	0.132E-02	0.117E-03	0.114E-01	0.852E-02 1589 32.53 25.94

Table 14 :

 14 Test 5 : Convergence tests for second order AB2 scheme, l 1 x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

	Rate	ρ -ρ ex l 1 x	p -p ex l 1 x	u -u ex l 1 x	ε -ε ex l 1 x
	AB1U	0.75	0.84	0.55	0.52
	AB1	1.07	1.16	0.58	0.72
	AB1M	1.23	1.33	0.84	0.91
	RK2U	0.81	0.88	0.64	0.58
	AB2U	0.75	0.82	0.57	0.54
	RK2	0.95	1.03	0.47	0.63
	AB2	0.97	1.04	0.47	0.65
	AB2M	1.29	1.37	1.01	0.95

Table 15 :

 15 Test 5 : Convergence tests for first and second order.

Table 16 :

 16 Test 6 : Convergence tests for first order AB1 scheme, l 1 x -norm at final time with respect to the averaged number of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

		1 x	p -p ex l 1 x	u -u ex l 1 x	ε -ε ex l 1 x	N T f	cpu g	cpu l
	1	0.645E+00 0.259E+02 0.104E+01 0.279E+03 200	1.30	1.30
	2	0.528E+00 0.192E+02 0.923E+00 0.208E+03 289	3.38	2.48
	3	0.409E+00 0.141E+02 0.801E+00 0.149E+03 435	9.16	6.21
	4	0.294E+00 0.104E+02 0.679E+00 0.102E+03 697 25.78 17.61
	5	0.197E+00 0.796E+01 0.588E+00 0.678E+02 1169 73.22 50.98
	L max	ρ -ρ ex l 1 x	p -p ex l 1 x	u -u ex l 1 x	ε -ε ex l 1 x	N T f	cpu g	cpu l
	1	0.548E+00 0.193E+02 0.872E+00 0.217E+03 200	1.31	1.36
	2	0.421E+00 0.132E+02 0.706E+00 0.156E+03 282	3.32	2.46
	3	0.301E+00 0.929E+01 0.566E+00 0.106E+03 418	8.82	5.64
	4	0.202E+00 0.672E+01 0.448E+00 0.700E+02 670 23.64 15.26
	5	0.127E+00 0.528E+01 0.396E+00 0.440E+02 1064 65.01 41.86

Table 17 :

 17 Test 6 : Convergence tests for second order AB2 scheme, l 1 x -norm at final time with respect to the averaged of cells N Lmax of ρ, u, p and ε, cpu g for global time stepping and cpu l local time stepping, N T f being the number of cells at the final time.

	Rate	ρ -ρ ex l 1 x	p -p ex l 1 x	u -u ex l 1 x	ε -ε ex l 1 x
	AB1U	0.35	0.60	0.41	0.47
	AB1	0.67	0.67	0.33	0.80
	AB1M	0.66	0.68	0.35	0.81
	RK2U	0.44	0.70	0.55	0.53
	AB2U	0.43	0.65	0.42	0.51
	RK2	0.90	0.75	0.48	0.98
	AB2	0.87	0.77	0.48	0.95
	AB2M	0.77	0.72	0.39	0.91

Table 18 :

 18 Test 6 : Convergence tests for first and second order.

Université de Toulon, IMATH, EA

2134, 83957 La Garde, France.December 13, 2013

Intel(R) Core(TM) i5-2500 CPU @ 3.30GHz

Acknowledgements

This work is supported by the ModTerCom project within the APEX program of the region Provence-Alpe-Côte d'Azur.