Gergö Barany
email: gergo.barany@inria.fr

Liveness-Driven Random Program Generation

Keywords: code generation, random testing, data-ow analysis, program optimization

Randomly generated programs are popular for testing compilers and program analysis tools, with hundreds of bugs in real-world C compilers found by random testing. However, existing random program generators may generate large amounts of dead code (computations whose result is never used). This leaves relatively little code to exercise a target compiler's more complex optimizations. To address this shortcoming, we introduce liveness-driven random program generation. In this approach the random program is constructed bottom-up, guided by a simultaneous structural data-ow analysis to ensure that the generator never generates dead code. The algorithm is implemented as a plugin for the Frama-C framework. We evaluate it in comparison to Csmith, the standard random C program generator. Our tool generates programs that compile to more machine code with a more complex instruction mix.

Motivation

Optimizing compilers for real-world programming languages are complex pieces of software. Compiler bugs may manifest in several ways: As compiler crashes, missed optimizations, or as silent miscompilations. The third category is especially serious as it may introduce bugs in correct programs. Such wrong-code bugs may invalidate all correctness guarantees provided by source-level verication of safety-critical (and other) software systems.

Two main avenues of work address these problems: compiler verication and compiler testing. Compiler verication has seen much research [START_REF] Dave | Compiler verication: A bibliography[END_REF], with Comp-Cert as a prominent example [START_REF] Leroy | Formal verication of a realistic compiler[END_REF]. However, such compilers have not entered the mainstream yet: Compiler verication is dicult and time-consuming, and veried compilers therefore perform fewer optimizations and target fewer CPU architectures than others.

A dierent approach is to test compilers in a way that instills condence. Standard compiler test suites exist for exercising C compilers, in particular for testing their conformance to various details of the standard [START_REF]ACVS ANSI/ISO/FIPS-160 C validation suite[END_REF][START_REF] Hall | The Plum Hall validation suite for C[END_REF]. In addition, randomized dierential testing has gained prominence in recent years. Compiling many random programs with various compilers and comparing the behaviors of the generated binaries can uncover input programs that cause compiler crashes or miscompilations. The best-known example of this approach is the work of Yang et al. on Csmith [START_REF] Yang | Finding and understanding bugs in C compilers[END_REF], a generator of random C programs. Csmith generates programs that are fully self-contained (including all their inputs in initialized global variables) and conform to the C standard by construction. If two compilers produce code that behaves dierently for a Csmith-generated program, one of the compilers must contain a miscompilation bug. Testing of C compilers with Csmith has uncovered hundreds of bugs in total, including crashes and miscompilations in every compiler under test. This included bugs in (unveried parts of) the CompCert veried C compiler [START_REF] Yang | Finding and understanding bugs in C compilers[END_REF].

This article describes a random generator of C programs developed for a project on nding missed optimizations in C compilers. Inspired by the successes of Csmith, in this project we generate random C programs, compile them using various compilers, then use custom tooling to search for possible optimizations in the resulting binaries. (The details are described in a separate paper [START_REF] Barany | Finding Missed Compiler Optimizations by Dierential Testing[END_REF].)

Starting with Csmith as our program generator, we found early on that it was not an optimal t for our intended use case: Despite generating realistic-looking programs with complex arithmetic expressions, accesses to global and local variables including through pointers, structures, and arrays, as well as nested loops and branches, it produces large amounts of dead code whose results are never used. (See our experiments in Section 4.) Dead code elimination, a standard part of every optimizing compiler, can thus remove large parts of the code generated by Csmith, leaving very little relevant code for the remaining more interesting optimizations. Csmith often generate functions of several hundred lines of code that are compiled to ten machine instructions, representing only a small fraction of the computations present on the source code level.

In this paper we address this problem with our new liveness-driven random generator ldrgen. Our tool uses liveness analysis during program generation to avoid generating dead code. In the following sections we describe liveness-driven program generation; the implementation of our tool as a Frama-C plugin; and its experimental evaluation, showing that ldrgen generates programs that compile to a larger amount of code and a more complex instruction mix than programs generated by Csmith.

Fully Live Programs

In this section we briey recall the basics of liveness analysis and then present our novel formulation as a set of structural inference rules.

Principles of Liveness Analysis

A variable is called live at a program point if the value it holds at that point may be read in the future, without an intervening redenition; otherwise, it is called dead. For example, in a code snippet like x = a + b; x = 0; return x;, the variable x is dead after the rst assignment but live after the second one (because it is used in the return statement). We can extend the notion of liveness from variables to the assignment statements dening them: An assignment v = e is live i the variable v is live just after this statement. The rst assignment to x above is dead, the second one is live. (Unfortunately, some authors use the term dead code to refer to unreachable code, as in if (false) x = y. These concepts are not the same; our use of the terms live and dead does not refer to reachability.)

Dead assignments without other side eects are useless and can be removed from the program. Even mildly optimizing compilers implement a dead code elimination pass that would completely remove the addition from the rst program fragment above. Our goal is to generate only live code, i. e., only code that does not contain any such opportunities for dead code elimination.

Liveness analysis is one of the classical data-ow analyses [START_REF] Nielson | Principles of Program Analysis[END_REF]. It is a backward, may analysis traditionally performed as backward xed-point iteration over a program's control-ow graph. A statement S in the control-ow graph has a live-in set S • and a live-out set S • which capture the sets of live variables before and after execution of S. Every statement S also has a transfer function f S relating these sets. Figure 1 shows the transfer functions for assignments, if statements, and while loops. Liveness information is noted on the edges of the control ow graph. In the equations, F V (e) denotes the set of all variables in expression e. The transfer function for an assignment v = e is said to kill the variable v and generate all the variables in e. The transfer functions for if statements and while loops kill nothing and generate the variables in their condition c. The live-in set S • of a loop S has a recursive dependency on itself: It is the same as the body's live-out set B • , on which the body's live-in set B • depends. The desired solution for S • is the least xed point of the system of equations, which can be found eciently by propagating data-ow information backwards around the loop.

v = e S • S • if (c) S • S • 1 S • 2 while (c) • • • S • S • B • B • S • S • = (S • \ {v}) ∪ F V (e) S • = S • 1 ∪ S • 2 ∪ F V (c) S • = S • ∪ B • ∪ F V (c)

Recognizing Fully Live programs

We will call a program fully live if all of the assignment statements it contains are live. This section develops an inference system characterizing fully live programs.

Figure 2 shows the abstract syntax of our programming language of interest, a subset of C function bodies without declarations. The language contains variables, constants, and all side-eect-free arithmetic and bitwise operators of C. Statements are assignments (the only source of side eects), return statements, if statements and general while loops. In contrast to C's concrete syntax, we view the semicolon ; as a statement separator, not a terminator. For now, the language does not include for loops, nor any structures, arrays, or pointers. All variables are considered local. Figure 3 shows a system of inference rules that characterize fully live programs. In these rules we use a notation similar to Hoare triples. A liveness triple

S • S S •
means that the variables in the set S • are live immediately before the statement S (live in), and the variables in S • are live immediately after S (live out). A program S is fully live i there is a set of variables S • such that the liveness triple S • S ∅ is derivable in the system. Intuitively, the system of inference rules encodes two things. First, the rules are an alternative presentation of the transfer functions of liveness analysis. A triple S • S S • that appears in a valid derivation corresponds to a dataow equation S • = f S (S • ) where f S is the transfer function for the statement S. For example, the transfer function f v=e for an assignment v = e is just f v=e (S • ) = (S • \ {v}) ∪ F V (e), as in the side condition of the Assign rule. Similarly, the Sequence rule encodes the composition of transfer functions, and the If rule encodes the split and join of data-ow information along dierent program paths.

Second, the other side conditions add a system of constraints to ensure full liveness: Any assignment statement appearing in a fully live program S (i. e., for which a derivation of S • S ∅ for some S • exists) is itself live. This follows directly from the Assign rule's side condition v ∈ S • . For example, a triple of the form

S • x = a; x = b S • Return {v} return v ∅ Skip S {} S Assign v ∈ S • S • = (S • \ {v}) ∪ F V (e) S • v = e S • Sequence S • 1 S1 S • 2 S • 2 S2 S • 2 S • 2 = ∅ S • 1 S1 ; S2 S • 2 If S • 1 S1 S • S • 2 S2 S • S • = S • 1 ∪ S • 2 ∪ F V (c) S1 = {} ∨ S2 = {} S • if (c) S1 else S2 S • While B • B B • B • = S • (minimal) S • = S • ∪ B • ∪ F V (c) S • = ∅ S • while (c) B S • Fig. 3.
System of inference rules for fully live programs can never be derived in the system because the rst of the two assignments is dead. The Sequence rule says that to derive this triple, there must be some intermediate set S of variables such that x ∈ S due to Assign on x = a while at the same time S = (S \ {x}) ∪ {b} due to Assign on x = b. This is a contradiction, so the derivation attempt must fail.

While the other rules are straight-forward, the While rule deserves some discussion. Unlike the two branches of the if statement, the whole loop's liveout set S • is not identical to the loop body's live-out set B • : Typically there are loop-carried dependences, i. e., cases where a variable is set on one iteration of the loop and its value is read on a later iteration. Such variables are live out of (and live into) the loop body, but if they are no longer used once the loop has terminated, they are not live out of the loop. When performing a derivation in the system, we must guess or calculate the set of these additional variables. Let f S denote the liveness transfer function corresponding to the loop body statement S. Then from the liveness triple

B • S B • we have B • = f S (B • ),
and the equation

B • = S • ∪ B • ∪ F V (c) means that B • is a xed point of the function λ B. (S • ∪ f S (B) ∪ F V (c)).
The minimality side condition additionally species that we are interested in the least xed point of this function. This least xed point exists and is unique [START_REF] Nielson | Principles of Program Analysis[END_REF].

In Figure 4 we illustrate the use of the inference system to prove full liveness of a program taking an input variable n (assumed to be non-negative) and returning the n-th Fibonacci number. We omit some details to focus on the analysis of the loop. Note that only the return variable a is live after the loop. However, the live-out set of the loop's body is {a, b, n}. This includes the return variable a and the variable n that is used in the loop condition. It also includes the variable b, which is the element computed by xed point iteration: The value of b at the end of the loop body will be used on the next loop iteration, if any. Conversely, if b were not live at some point in the loop body, our inference system would not allow derivation of a triple for the assignment b = t. Indeed, all assignments in the loop body satisfy the condition that they dene variables that are live after the assignment. That is, this program is fully live by our denition.

{a, b, n} t = a + b {b, n, t} {b, n, t} a = b {a, n, t} {a, n, t} b = t {a, b, n} {a, b, n} n = n -1 {a, b, n} • • • • • • • • • {a, b, n} t = a + b; a = b; b = t; n = n -1 {a, b, n} {a, b, n} while (n > 0) { t = a + b; a = b; b = t; n = n -1 } {a} . . . {n} a = 0; b = 1; while (n > 0) { t = a + b; a = b; b = t; n = n -1 }; return a ∅

Limitations of the System

Note that fully live programs may still contain opportunities for simple optimizations that remove code that does not have interesting eects. For example, programs accepted by the inference system above may contain fragments like if (0) { ... } else { ... } where one of the branches of the if statement is unreachable and thus irrelevant; or assignments like y = ...; x = y -y; where the computation for the value of y is irrelevant for x's nal value of 0. Our inference rules do not consider the semantics of the code in enough detail to exclude such cases.

Our claims with regards to full liveness are relative to a purely syntactic notion of liveness that does not consider such semantic issues. In particular, we cannot guarantee that the liveness analysis embedded in these rules is equivalent to liveness analysis as performed by any given compiler. Any other analyses or transformations performed by the compiler before liveness analysis may inuence the results, typically making the compiler's results more precise than ours.

As our experimental results in Section 4 show, our generator performs well nonetheless, so we can leave renements of the system for future work.

Generating Fully Live Programs

The inference rules can be translated almost directly into an executable random (or exhaustive) generator of fully live programs. Like traditional liveness analysis, generation proceeds backwards, i. e., in the direction opposite control ow.

The side conditions of the inference rules ensure that a fully live program always ends in a return statement, as no other statement may have an empty live-out set. The generator can thus start by picking a random program variable v and generating a statement return v with live-in set L = {v}. It then iteratively prepends random statements S to the current program fragment and updates the live-in set according to f S (L). The possibilities for the generation of S are guided by L. In particular, if the generator decides to generate an assignment statement, the target variable v must be an element of L at that point. Conversely, if L ever becomes empty, generation of the current block of code must stop at that point: Any code preceding that point would be dead. Figure 5 shows pseudocode of such a generator in an OCaml-like functional language.

Every statement generation function takes a live variable set L (representing the live-out set S • of the statement S to be generated) and returns a pair of a newly generated statement and an updated live variable set according to the statement's transfer function. We iterate statement generation and collect a list of statements forming a block. In this presentation we omit helper functions such as the ones for generating random variables and expressions.

As before, the handling of loops merits more discussion. Just as the inference system needs a minimal set B • containing new variables that are live into and out of the loop body, the random generator must synthesize such a set of variables. But here the problem is more dicult: During inference, we can start with an initial live-out set and derive the eventual live-in set by xed-point iteration. During code generation this is not possible since we cannot analyze the loop body before we have constructed it. Instead, we rst generate a random set of newly live variables and let this choice guide generation of the loop's code.

This solution relies on the following observation: The new variables we are interested in are ones that are dened before the loop, may be dened on some loop iteration, and then used on some later iteration. In the example of Figure 4, this is the case for variable b, which stores the next Fibonacci number for assignment to a on the subsequent iteration. Using the names in the While rule of Figure 3, the set of these `new' variables is

B = B • \ (S • ∪ F V (c)). It follows that B ⊆ B • , i. e., every b ∈ B is live into the loop.
To generate a fully live loop, we choose a random set B of new variables and generate a loop body in a way that ensures that

B • = S • ∪ B ∪ F V (c)
is a least xed point of the loop. For this we add B to the live variable set L before generating a loop body block (along with the live variables F V (c) generated by the loop condition). The generated loop body may or may not dene or use variables in B , and thus these variables may or may not be live into the generated block. However, we must force them to be used in the loop body and be live into the body: If there were some b in the generated B that is not live into the loop body, we would violate the condition On the other hand, if b were live into the loop body but did not have a use anywhere in the body, then B • would not be minimal and hence B • would not be a least xed point of the constraint system.

("v = e", (L \ {v}) ∪ F V (e)) let branch L = let (t, L1) = random_statements L []; let (f, L2) = random_statements L []; let c = random_expression (); ("if (c) t else f ", L1 ∪ L2 ∪ F V (c))
To ensure a correct, minimal solution, we therefore nd the set V of all b B that are not in L or that are in L but have no use in the generated loop body. We pick a random live variable v ∈ L and prepend an assignment v = e to the generated loop body, where e is an expression containing all the variables in V . This nal loop body ensures that all variables in B are live into it and used in it, hence ensuring that

B • = S • ∪ B ∪ F V (c
) is a least xed point of the loop's liveness constraint system.

A small detail not illustrated in the pseudocode is the case when L is empty at the beginning of the generated loop body. This can only be the case if the rst statement in the body is an assignment of a constant expression (i. e., not using any variables) to some variable v, since such assignments are the only statements that can remove variables from the live variable set without adding any new ones. In this case, we replace this assignment's right-hand side with the expression e generated as above.

Implementation

We have implemented the random program generation algorithm sketched above as a plugin for Frama-C [START_REF] Kirchner | Frama-C: A software analysis perspective[END_REF]. Frama-C is a general, extensible framework for source-level analysis and transformation of C programs. It is written in OCaml and can be extended with plugins written in that language. For this work, we do not need any advanced Frama-C features but benet from its AST (abstract syntax tree) type denitions, utilities for managing variables and constructing AST fragments, and its pretty-printer for outputting the generated AST as C source code. As these general parts are provided by Frama-C, ldrgen itself can be quite small: It consists of only about 600 lines of generator code, plus some utilities and conguration.

ldrgen is free software, available at https://github.com/gergo-/ldrgen.

Random Generation

The core of the generator has the same structure as the pseudocode in Figure 5.

After generating an empty function denition and a return statement, it lls in the function's body by generating a fully live sequence of statements as in the pseudocode. Statements are represented by AST fragments; we never need to worry about generating actual C syntax. The current version of ldrgen always generates a single function. Random expressions are generated by choosing an operator among the arithmetic operators available in C and recursively generating the appropriate number of operand expressions. At this point, C's type system becomes relevant; if needed, we insert type casts to ensure all operands of an operator have the same type. Type casts are also needed in some other cases: Bitwise operators and the modulo operator cannot be applied to oating-point numbers in C, so we insert conversions to integer types in such cases.

Many C operators may invoke undened behavior when applied to inappropriate values. Two examples are division by zero and signed integer arithmetic overow. Unlike Csmith in its default mode, ldrgen does not try to guard against such undened except for two cases that compilers have repeatedly warned us about: We clamp the right-hand-side operands of bit-shift operations to the bit size of the expression on the left-hand-side, and we always generate division and modulo operations of the form e 1 /(e 2 +c) for some constant c instead of just e 1 /e 2 . The idea behind this is that e 2 +c is less likely to evaluate to zero than a random expression in general. This approach is primitive, but we have found it to work well in practice.

Leaves of expressions are constants or variables. For constant literals we simply generate a random number. For a variable use we either use a previously used variable or generate a new one. Variables generated in this way may be local variables or function parameters. Both can be used in expressions, but we only generate assignments to locals, not to parameters. Some of the generator's choices are weighted by manually chosen parameters to ensure generation of somewhat more realistic-looking programs. For example, we prefer generation of basic arithmetic operations to bitwise operators. We also ensure that loop and branch conditions are not constant expressions, i. e., that they contain at least one variable. In order to avoid trivial non-termination issues, we also ensure that every loop body's nal statement is an assignment to some variable that occurs in the loop condition. If there were no modication of any of these variables at all, a loop once entered could never terminate. Even so, termination is not at all guaranteed.

Bottom-up generation of the function's body may stop if there are no more live variables, or if a user-dened limit is reached. In this latter case, there may remain live local variables at the start of the function's body. Their liveness means that they may be used without being assigned to, so we must ensure that they are initialized. We therefore nalize the function denition by initializing all such live-in variables to constants or to the values of function parameters.

Conguration

The generator's behavior may be tuned using command-line arguments. These may specify features of the sub-language of C that is used. For example, the user may request the generation of code that only uses integer types, or only oating-point types. They may also specify that no bitwise operations or no divisions should be generated, and whether loops may be generated. Other ags specify structural properties: The maximal number of statements per block, and the maximal nesting depths of statements and expressions.

ldrgen's random generation uses OCaml's standard pseudorandom number generator, which can be seeded with a random seed or with a seed value specied as a command line argument. Invoking a given version of ldrgen with a xed set of arguments and a xed seed thus always gives the same reproducible result.

Extensions to the Basic Model

We describe two extensions to the core language of Figure 2 that are already implemented in ldrgen: very limited uses of pointers and for loops over arrays.

First, in addition to the arithmetic types used so far, we can generate function parameters of type T * (pointer to T) for some arithmetic type A parameter p of such type can be used in generated code as *p. We currently do not generate assignments to such dereferenced pointers, nor any pointer arithmetic.

Second, we want to generate arrays and restricted forms of loops over them in order to exercise loop optimizations such as unrolling or vectorization. For this we generate pointer arguments T *arr which are only used in for loops of the following form:

v = . . .; for (unsigned int i = 0; i < N; i++) { v = v • f (arr[i]); }
Here, N is a global variable considered to hold the array's size, f (arr[i]) is a random expression involving arr[i], and • is a randomly chosen binary arithmetic operator. This loop pattern implements a map-reduce operation, mapping some function f over the array and reducing (folding) the result with •. It is currently the only kind of for loop we implement, but this would be easy to generalize.

Similar forms of loops are already generated by Csmith, but there their results are virtually never used. Csmith's loops are therefore completely eliminated by compilers instead of being vectorized. In ldrgen, we choose a loop result variable v that is live after the loop to ensure that it is used, and the loop exposed to the compiler's loop optimizer.

Future Extensions

In the future, we are planning to extend ldrgen to generate structure types and allow the use of their members.

In the longer term, ldrgen will also be extended to support programs consisting of several random functions which may call each other. We are not planning to support non-structured control ow using goto. The more structured break and continue statements might eventually be supported, but this is not a priority as they complicate the structural liveness analysis.

Evaluation

The design goal of ldrgen was to have a random program generator that exposes as much interesting code as possible to all passes of the compiler under test; recall that we found Csmith-generated code to contain much dead code which is never seen by many parts of the compiler because it can be optimized away early on.

We will have achieved our goal if, for comparable amounts of generated C code, ldrgen's output results in more, and ideally more varied, assembly code than Csmith's output. We therefore compare the two generators along these lines. We do not claim superiority to Csmith in any other regard, especially not concerning its power to nd subtle miscompilation bugs. Csmith covers a larger subset of the C language than ldrgen and, in its default mode, carefully ensures that its output is well-dened according to the C standard. ldrgen tries to guard against oversized shifts and divisions by zero (see Section 3.1) but does not otherwise guarantee well-denedness.

Csmith is designed to run complete, self-contained applications consisting of several functions, driven by a main function. In contrast, ldrgen only generates individual functions without a driver. However, Csmith's many conguration options allow us to ask it to generate les consisting only a single function without main. 1Table 1 presents our experimental results for 1000 programs each generated by Csmith and ldrgen. We investigate three characteristics of the generated programs: lines of C code, number of instructions in the generated code, and number of unique opcodes in the generated code. In all cases, the C code was compiled to x86-64 machine code using GCC 5.4.0 with optimization setting -O3. For each characteristic, the table shows the total over the 1000 les as well as the minimum, median, and maximum values. (In cases where the median is not unique, we chose the arithmetic mean of the two closest values.) Our command line ags for Csmith were chosen in order to generate comparable numbers of lines of C code to ldrgen. In fact Table 1 shows that it generates somewhat more, but these numbers are dicult to compare precisely because Csmith-generated code tends to contain many initializers for global variables; ldrgen does not generate any global variables at all. We believe that the settings we chose allow a fair comparison of the generators.

Next we compare the number of instructions (executable code only, excluding static data, assembler directives etc.) emitted by the compiler for the generated source les. ldrgen was designed to increase this number compared to Csmith, and the table shows that we have succeeded: While on average Csmith's output compiles to a single machine instruction per ten lines of code, ldrgen's output has almost three instructions per single line of source code. Overall, ldrgen-generated programs compile to about 20 times as much machine code as Csmith-generated programs of comparable size. We can also see that the distribution for Csmith is highly skewed: The median shows that at least half of the functions generated by Csmith compile to 15 instructions or fewer. This also conrms our initial, more informal observation that Csmith-generated code tends to contain large amounts of dead code. ldrgen manages to generate code with a less skewed distribution, and in particular with generally higher numbers of emitted instructions.

On a side note, we remark that both Csmith and ldrgen sometimes generate functions that compile to a single machine instruction. Inspection showed that this happens in cases where the compiler recognizes that a function ends up in an innite loop without externally visible side eects. Such functions are then compiled into a single unconditional jump instruction looping back to itself. Many other functions compile to two instructions, typically some simple operation on a function argument or a constant followed by a return. It would be dicult to completely avoid generating innite loops, but comparatively easy (at least within ldrgen) to avoid generating functions that return after a single operation. For both Csmith and ldrgen, about 10 % of all cases fall into one of these trivial categories (with Csmith producing fewer innite loops).

We analyze the coverage of the instruction set in the generated code by looking at the number of dierent opcodes generated. Here, too, we see that individual functions generated by ldrgen have a more varied instruction mix than functions generated by Csmith: Even the median for ldrgen is higher than the maximum for Csmith. Totaling over all the machine code in 1000 functions, we see that Csmith-generated code compiles to a mix of 146 dierent opcodes, while ldrgen-generated code contains 204 dierent opcodes, an increase in instruction set coverage of 40 %. Inspection of the sets of opcodes shows that this dierence is almost entirely due to various vector (SIMD) arithmetic instructions generated for ldrgen's code. Compiling to such instructions was the goal of adding for loops over arrays to ldrgen. Manual inspection of some cases shows that such loops are indeed the origin of these instructions. As noted above, such loops are also generated by Csmith, but their results are almost never used, so they do not appear in the compiled code. Disabling generation of for loops in ldrgen brings its total number of unique instructions down to 147, comparable to Csmith.

One of the few opcodes emitted for Csmith-generated code but not for ldrgen are call instructions to memcpy which are sometimes generated by compilers for structure copies. ldrgen currently does not generate structures at all.

Finally, in Table 2 we compare the speed of the two generators. Generating the 1000 les each analyzed above took 871 seconds with Csmith and 124 seconds with ldrgen (Csmith backtracks if it nds that it has generated unsafe code). Csmith generates about 527 lines of C code per second, with ldrgen generating 3140 (about 6 × more). With respect to nal machine code, Csmith-generated code compiles to about 52 instructions per second of generation time, whereas ldrgen produces 8563 (about 160 ×). These numbers do not include the time taken by the compiler; compiling all 1000 les for each generator takes about 46.6 seconds for Csmith and 80.3 seconds for ldrgen. The best-known random program generator is Csmith [START_REF] Yang | Finding and understanding bugs in C compilers[END_REF], based on an earlier system called randprog [START_REF] Eide | Volatiles are miscompiled, and what to do about it[END_REF][START_REF] Turner | Random C program generator[END_REF]. Csmith generates complete, self-contained programs that take all their input from initialized global variables and compute an output consisting of a hash over the values of all global variables at the end of execution. The generator is designed to only generate programs with welldened semantics: Operations that may be undened in C, such as overowing signed integer arithmetic, are guarded by conditionals that exclude undened cases (these guards can be disabled, and we disabled them for the experiments reported above). Like ldrgen, Csmith performs data-ow analysis during generation, although the details dier due to the diering design goals. Csmith's forward analysis computes points-to facts and uses them for safety checks. If the checks fail, Csmith backtracks, deleting code it generated until a safe state is reached again. In contrast, ldrgen's data-ow analysis only deals with liveness, and ldrgen never backtracks: Full liveness of variables in loops is ensured by construction. Csmith generates a larger subset of C than current or currently planned versions of ldrgen, including unstructured control ow and less restricted use of pointers. Due to Csmith's forward generation and data-ow analysis, it does not appear possible to directly integrate our backward liveness-driven approach in Csmith. A `best eort' approach that would not give full liveness guarantees might extend Csmith with a forwards reaching denitions analysis [START_REF] Nielson | Principles of Program Analysis[END_REF]. At each program point, when generating an expression, Csmith could then prefer to use previously dened variables over fresh variables. Alternatively, it might be extended with a mode that computes a function's return value from all otherwise unused variables in the program, articially making them live.

Csmith has been used to nd hundreds of bugs in C compilers when compiling the programs it generates [START_REF] Yang | Finding and understanding bugs in C compilers[END_REF]. It has also been used as the basis of mutationbased systems, where Csmith's output was modied using other tools to provoke compiler bugs [START_REF] Le | Compiler validation via equivalence modulo inputs[END_REF]. The CLsmith tool derived from Csmith has been used to nd many bugs in OpenCL compilers [START_REF] Lidbury | Many-core compiler fuzzing[END_REF].

Other Random Generators

Other notable random generators include Orange3 [START_REF] Nagai | Reinforcing random testing of arithmetic optimization of C compilers by scaling up size and number of expressions[END_REF] for C, jsfunfuzz [START_REF] Ruderman | [END_REF] for JavaScript, and etester [START_REF] Midtgaard | Eect-driven quickchecking of compilers[END_REF] for OCaml. Our generator following a system of inference rules is similar to a generator based on Haskell's typing rules [START_REF] Paªka | Testing an optimising compiler by generating random lambda terms[END_REF].

The JTT program generator [START_REF] Zhao | Automated test program generation for an industrial optimizing compiler[END_REF] is aimed at testing compiler optimizations. It uses a model-based approach, with generation guided by test scripts. These scripts contain code templates and temporal logic specications of the optimizations to be tested. For example, the authors specify opportunities for dead code elimination as cases where a variable is assigned, then assigned again before being used. The test script contains a temporal logic formula expressing this pattern and the test condition that the compiler should eliminate the rst assignment. Using this script, JTT generates test programs containing this pattern.

Other work specically aimed at testing and comparing program verication tools generates code from randomly generated LTL formulae [START_REF] Steen | Property-driven benchmark generation: synthesizing programs of realistic structure[END_REF]. The generated code is guaranteed to satisfy the specied temporal properties.

Structural Data-Flow Analysis

Our formulation of liveness analysis as set of structural inference rules is inspired by formulations of interval-based data-ow analysis where reducible programs are decomposed into components called intervals, and analysis data is eciently propagated among the intervals [START_REF] Allen | Control ow analysis[END_REF][START_REF] Cocke | Global common subexpression elimination[END_REF][START_REF] Graham | A fast and usually linear algorithm for global ow analysis[END_REF]. [START_REF] Graham | A fast and usually linear algorithm for global ow analysis[END_REF]

Conclusions

We presented ldrgen, a new generator of random C programs designed for testing C compilers. In contrast to Csmith, the dominant player in this eld, ldrgen is driven by liveness analysis to avoid generating dead code. We designed an inference system to capture our liveness analysis and implemented its rules as an excutable program generation system.

ldrgen is implemented as a plugin for the Frama-C framework. Our evaluation of ldrgen in comparison to Csmith shows that we have achieved our goal of generating C code that compiles to larger amounts of machine code with a more interesting instruction mix, including many SIMD instructions. We are actively using ldrgen in a project on nding missed optimizations in compilers [START_REF] Barany | Finding Missed Compiler Optimizations by Dierential Testing[END_REF]. Because it is able to exercise loop optimizations not usually addressed by Csmith, it may also be useful for nding correctness bugs in these optimizations.

Fig. 1 .

 1 Fig. 1. Liveness transfer functions

vFig. 2 .

 2 Fig. 2. Abstract syntax of a C-like programming language

Fig. 4 .

 4 Fig. 4. Example derivation proving full liveness.

  let loop L = ( * See main text for explanation. * ) let c = random_expression (); let B = random_variable_set (); let (code, L') = random_statements (L ∪ B ∪ F V (c)) []; let V = {b ∈ B | b / ∈ L or b not used in S}; if V = ∅ then let e = random_expression_on_variables V; let v = random_select L ; let code' = "v = e" :: code; ("while (c) code ", (L \ {v}) ∪ V ∪ else ("while (c) code", L ∪ L) ( * Start generation with the terminating return statement. * ) let v = random_variable () let (code, L) = random_statements {v} ["return v"]

Fig. 5 .

 5 Fig. 5. Pseudocode of a liveness-driven random program generator.

  B ⊆ B • established above.

	let random_statements L code =
	if L = ∅ then
	let (S, L ) = random_statement L;
	random_statements L (S :: code)
	else
	(code, L)
	let random_statement L =
	let statement_generator = random_select [assignment; branch; loop];
	statement_generator L
	let assignment L =
	let v = random_select L;

let e = random_expression ();

Table 1 .

 1 Comparison of code generated by Csmith and ldrgen in 1000 runs each.

		generator min median max	total
	lines of code	Csmith	25	368.5 2953	459021
		ldrgen	12	411.5 1003 389939
	instructions	Csmith	1	15.0 1006	45606
		ldrgen	1	952.5 4420 1063503
	unique opcodes Csmith	1	8	74	146
		ldrgen	1	95 124	204

Table 2 .

 2 Comparison of the time to generate 1000 les.

	generator time (sec) lines/sec instrs/sec
	Csmith	871	527	52.4
	ldrgen	124	3140	8562.8
	5 Related Work
	5.1 Csmith

The concrete ags we used were nomain float max-funcs 1 no-safe-math max-block-size 8 concise.

Acknowledgments. The author would like to thank the anonymous reviewers, John Regehr, and Gabriel Scherer for insightful comments on earlier versions of this paper. This research was partially supported by ITEA 3 project no. 14014, ASSUME.


Previously in

Until recently, lipid droplets (LDs) were just considered as neutral lipid storage sites, thus providing energy through lipolysis and β-oxidation pathway when required during stressful conditions [START_REF] Farese | Lipid droplets finally get a little R-E-S-P-E-C-T[END_REF].

Nowadays, it's clearly established that LDs are well-organized and extremely dynamic cellular organelles, conserved in eubacteria, fungi, plants and animals, where they are essential for lipid homeostasis and energy maintenance [START_REF] Murphy | The biogenesis and functions of lipid bodies in animals, plants and microorganisms[END_REF]. Alteration of these two critical physiological processes can lead to important metabolic disorders [START_REF] Welte | Expanding roles for lipid droplets[END_REF]. The "conventional snapshot representation" of LDs remains a suitable model for defining their general composition which is based on a central organic core of neutral lipids (mainly triacylglycerol (TAG) and sterol esters) surrounded by a monolayer of phospholipids [START_REF] Tauchi-Sato | The surface of lipid droplets is a phospholipid monolayer with a unique Fatty Acid composition[END_REF], associated with a wide range of structural [START_REF] Greenberg | Perilipin, a major hormonally regulated adipocyte-specific phosphoprotein associated with the periphery of lipid storage droplets[END_REF], enzymatic [START_REF] Sztalryd | Perilipin A is essential for the translocation of hormone-sensitive lipase during lipolytic activation[END_REF] and membranetrafficking proteins [START_REF] Liu | Chinese hamster ovary K2 cell lipid droplets appear to be metabolic organelles involved in membrane traffic[END_REF]. However, we also know that proteins and lipid species of such structures can be extremely diversified depending on the different cell types or metabolic status. Interestingly, LDs anabolism and catabolism are well-balanced and tightly controlled biological mechanisms involving a large number of actors at both transcriptional, translational and post-translational levels [START_REF] Lass | Lipolysis -a highly regulated multi-enzyme complex mediates the catabolism of cellular fat stores[END_REF]. TAGcontaining LDs breakdown is achieved during extended starvation period or enhanced energy demand, and this phenomenon is mainly mediated by three distinct lipolytic enzymes (i.e. the adipose triglyceride lipase (ATGL), the hormone-sensitive lipase (HSL) and the monoglyceride lipase (MGL)) that act sequentially to finally generate free fatty acids (FFA) and glycerol molecules [START_REF] Lass | Lipolysis -a highly regulated multi-enzyme complex mediates the catabolism of cellular fat stores[END_REF]. Since ATGL is catalyzing the first step of this essential lipolytic pathway, it's crucial to fully define physiological function(s) and structural properties of this protein but also to obtain further insights onto the regulatory mechanisms governing ATGL action towards LDs [START_REF] Zechner | Adipose triglyceride lipase and the lipolytic catabolism of cellular fat stores[END_REF].

In a fascinating way, ATGL was discovered fourteen years ago by three independent groups at the same time [START_REF] Zimmermann | Fat mobilization in adipose tissue is promoted by adipose triglyceride lipase[END_REF][START_REF] Jenkins | Identification, cloning, expression, and purification of three novel human calcium-independent phospholipase A2 family members possessing triacylglycerol lipase and acylglycerol transacylase activities[END_REF][START_REF] Villena | Desnutrin, an adipocyte gene encoding a novel patatin domaincontaining protein, is induced by fasting and glucocorticoids: ectopic expression of desnutrin increases triglyceride hydrolysis[END_REF]. Reports described that the protein was displaying a strong TAG-hydrolase activity in both in-vitro and ex-vivo experimental conditions [START_REF] Zimmermann | Fat mobilization in adipose tissue is promoted by adipose triglyceride lipase[END_REF][START_REF] Jenkins | Identification, cloning, expression, and purification of three novel human calcium-independent phospholipase A2 family members possessing triacylglycerol lipase and acylglycerol transacylase activities[END_REF][START_REF] Villena | Desnutrin, an adipocyte gene encoding a novel patatin domaincontaining protein, is induced by fasting and glucocorticoids: ectopic expression of desnutrin increases triglyceride hydrolysis[END_REF]. Moreover, Atgl gene was highly expressed within adipose tissues, and to a lesser level in liver, spleen, kidney, heart and skeletal muscle [START_REF] Zimmermann | Fat mobilization in adipose tissue is promoted by adipose triglyceride lipase[END_REF][START_REF] Villena | Desnutrin, an adipocyte gene encoding a novel patatin domaincontaining protein, is induced by fasting and glucocorticoids: ectopic expression of desnutrin increases triglyceride hydrolysis[END_REF]. To better understand the physiological role of this 54 kDa protein in lipid homeostasis, an Atgl -/-mutant mouse was generated, and study of this mouse permitted to obtained the first evidences that the ATGL protein is playing an essential role in TAG hydrolysis in-vivo [START_REF] Haemmerle | Defective lipolysis and altered energy metabolism in mice lacking adipose triglyceride lipase[END_REF].

Indeed, several phenotypes related to lipid metabolism disorder were easily observed within an Atgl null mutant such as increase in body weight, fat mass, fat accumulation in non-adipose tissues, and also a greater resistance to glucose and insulin. In addition, this deficiency rapidly triggered TAG accumulation within cardiac muscle thus leading to cardiac dysfunction and premature death [START_REF] Zechner | Adipose triglyceride lipase and the lipolytic catabolism of cellular fat stores[END_REF][START_REF] Haemmerle | Defective lipolysis and altered energy metabolism in mice lacking adipose triglyceride lipase[END_REF].

All these findings suggested that a new essential component was involved in central lipid metabolism, and thus opening new perspectives to better control lipid metabolic disorder in patients. Among the proteins involved in this regulation process, Plin1 (one of the five members of the perilipin family (Plin1-Plin5)) CGI-58, G0S2 and PEDF are probably the most important factors impacting ATGL activity [START_REF] Cerk | Adipose Triglyceride Lipase Regulation: An Overview[END_REF]. Plin1 is mastering the switch from basal to stimulated lipolysis, and this is mediated by its C-terminal domain which sequesters the CGI-58 protein and so prevents the action of the ATGL [START_REF] Lass | Adipose triglyceride lipase-mediated lipolysis of cellular fat stores is activated by CGI-58 and defective in Chanarin-Dorfman Syndrome[END_REF]. Upon specific hormonal stimulation Plin1 is phosphorylated and releases the CGI-58 co-activator protein which binds to ATGL and leads to its translocation at the LDs surface [START_REF] Lass | Adipose triglyceride lipase-mediated lipolysis of cellular fat stores is activated by CGI-58 and defective in Chanarin-Dorfman Syndrome[END_REF]. Point mutations, insertions or deletions within the cgi-58 gene trigger a drastic neutral lipid storage disorder also called Chanarin-Dorfman syndrome thus emphasizing the role of CGI-58 protein in LDs degradation [START_REF] Lass | Adipose triglyceride lipase-mediated lipolysis of cellular fat stores is activated by CGI-58 and defective in Chanarin-Dorfman Syndrome[END_REF][START_REF] Lefevre | Mutations in CGI-58, the gene encoding a new protein of the esterase/lipase/thioesterase subfamily, in Chanarin-Dorfman syndrome[END_REF]. In addition to CGI-58, the PEDF protein is also known for interacting with ATGL and stimulating lipolysis in adipocytes [START_REF] Borg | Pigment epithelium-derived factor regulates lipid metabolism via adipose triglyceride lipase[END_REF]. In contrast to CGI-58 and PEDF which are activators, G0S2 protein is negatively regulating ATGL activity in both in-vitro and in-vivo conditions where overexpression leads to an almost identical phenotype than an Atgl -/-mutant [START_REF] Yang | The G(0)/G(1) switch gene 2 regulates adipose lipolysis through association with adipose triglyceride lipase[END_REF][START_REF] Heckmann | Defective adipose lipolysis and altered global energy metabolism in mice with adipose overexpression of the lipolytic inhibitor G0/G1 switch gene 2 (G0S2)[END_REF][START_REF] Cerk | A peptide derived from G0/G1 switch gene 2 acts as noncompetitive inhibitor of adipose triglyceride lipase[END_REF].

Another important part of their manuscript was dedicated to the inhibition of ATGL activity by either natural or synthetic small molecules [START_REF] Cerk | Adipose Triglyceride Lipase Regulation: An Overview[END_REF]. One of the main inhibitory mechanism towards ATGL and also HSL is mediated by acyl-CoA availability within the cells [START_REF] Severson | Inhibition of the hormone-sensitive lipase in adipose tissue by long-chain fatty acyl coenzyme A[END_REF][START_REF] Nagy | Adipose triglyceride lipase activity is inhibited by long-chain acylcoenzyme A[END_REF]. Indeed, an increase level of such molecules drastically impairs LDs-associated lipases activities and could directly contribute to the feedback inhibition of lipolysis. Finally, during high throughput screening of chemical compounds, one specific synthetic molecule has been identified as powerful inhibitor of ATGL. This compound named Atglistatin selectively inhibits the mouse ATGL activity at a micro/nanomolar range and drastically reduce TAG and FFA plasma level [START_REF] Mayer | Development of small-molecule inhibitors targeting adipose triglyceride lipase[END_REF].

A large number of open questions needs to be further investigated regarding these dynamic interactions and new approaches a currently developed to better understand such mechanisms. For example, by generating translational fusions between the APEX2 protein and either the Plin2 or the ATGL protein, Bersuker et al., recently define with an high confidence a dynamic LD proteome in human cells [START_REF] Bersuker | A Proximity Labeling Strategy Provides Insights into the Composition and Dynamics of Lipid Droplet Proteomes[END_REF]. By using this powerful proximity labelling strategy, they were able to identify new LDs-associated proteins but also to describe new potential interactions between structural, enzymatic and membrane-trafficking proteins [START_REF] Bersuker | A Proximity Labeling Strategy Provides Insights into the Composition and Dynamics of Lipid Droplet Proteomes[END_REF]. LDs metabolism plays an important role in several diseases, such as obesity, atherosclerosis, metabolic syndrome, neurodegenerative diseases and mitochondrial disorders, which often lead to diabetes and cardiovascular complications. Altogether, these information summarized by Cerk et al., demonstrated that understanding the molecular mechanism of ATGL action and its regulation are crucial to further developed potent molecules for the treatment of such neutral lipid storage disorder.

  Current Protein & Peptide Science, Cerk et al. (2018) presented an updated review of concepts and knowledge regarding structure, function and regulatory mechanism of the adipose triglyceride lipase (ATGL), one of the key-enzyme required for intracellular lipolysis.

  These information prompted several teams to further investigate, during the last decade, the potential regulatory mechanisms involved in ATGL activity and to date, a large number of ATGL posttranslational modifications have been identified. In this context Cerk, Wechselberger and Oberer report in their recent paper published in Current Protein & Peptide Science, 2018;19(2):221-233, an updated and nice overview of knowledge regarding ATGL function with special focus onto these post-translational regulatory mechanisms impacting its TAG-hydrolase activity during LDs breakdown [14].

Financial Supports

PS received financial support for his PhD fellowship from the Ministère Français de l'Enseignement Supérieur, de la Recherche et de l'Innovation


Filippo Masi
email: filippo.masi@enpc.fr

Paolo Maria Mariano
email: paolo.mariano@unifi.it

Paolo Vannucci
email: paolo.vannucci@uvsq.it

Blast actions in aircrafts: an integrated methodology for designing protection devices

Keywords: Continuum mechanics, blast actions, shock waves, explosions, aircrafts

We introduce a numerical approach for evaluating the behavior of aircraft fuselages subjected to internal explosions. At variance with available approaches, we consider fluid-structure interactions through a novel integrated methodology able to take into account both stress at cruise altitude and blast fast-dynamics in interaction with pressurization not considered as a static load. The precision of the proposed numerical procedure allows us to foresee a nonstandard composite protective device.

Introduction

We record non-accidental airplane internal explosions since 1933, in the cargo hold of a Boeing 247D; the explosive was nitroglycerin. Table 1 shows a list of similar subsequent events. Until the explosion induced in the Boeing 707-124 in 1962, explosive devices placed in the baggage compartment plaid role. Then, in-cabin explosions became dominant. First, in 1987 liquid explosives were used and rapidly replaced by plastic explosives placed inside shoes, laptops, and other devices. Prevention based on a screening before boarding contrasted so far these actions. Gillen and Morrison [START_REF] Gillen | Aviation security: Costing, pricing, finance and performance[END_REF] report a comparative study of European total expenditures on aviation security: 5.7 billion euros in 2011.

In this context, the idea of the so-called unit load devices (ULD) emerged. It is a design of luggage containers with the aim of absorbing energy from an in-cargo explosion. Examples are ULDs made of fiber-reinforced composites [START_REF] Mlakar | Hardened Unit Load Devices[END_REF] and bilayer hardened luggage containers. In the latter case, the inner layer of lightweight foam captures debris, the outer layer mitigates pressure [START_REF] Sanai | Hardened Luggage Container[END_REF]. Usual protections (see e.g. [START_REF] Fray | Venting a space to relieve pressure generated by an explosion[END_REF] and [START_REF] Hararat-Tehrani | Decompression panel for a separation device in an aircraft[END_REF]) consist of blowout panels designed to be weaker than the sorrounding airframe. During an in-cabin explosion, blowout panels fail with consequent pressure decrement and possibly controlled fuselage failure. At a cruise altitude, pressurization and gravity play a non-negligible role together with the inertia of the rigid-body component of the airplane motion.

Standard experiments on the overall mechanical behavior of fuselages usually deal with a fatigue design. Those involving explosions commonly exploit an aircraft at ground, loaded just by gravity (see [START_REF] Settles | Shock waves in aviation security and safety[END_REF]). The experiment described in reference [START_REF] Backman | Blast Response of a Pressurized Aircraft Fuselage Measured Using High-Speed Image Correlation[END_REF] considers a partial pressurization in a Boeing 727, while those reported in reference [START_REF] Veldam | Response of pre-pressurized reinforced plates under blast loading[END_REF] refer just to a plane panel with a preceding pressurization. Large-scale effects afflict fuselage dynamics.

We record [START_REF] Soutis | Modelling the structural response of GLARE panels to blast load[END_REF] attempts of designing reinforced plates made of Aluminum-based alloys or glass-reinforced Aluminum (GLARE) by taking into account blast actions. A question not yet largely investigated is, however, the behavior of the aircraft in its whole.

There are computational analyses of blast actions on fuselages, based on a Coupled-Eulerian-Lagrangian (CEL) approach. In particular A. Dacks and J. Toczyski [START_REF] Dacko | Vulnerability analysis of aircraft fuselage subjected to internal explosion[END_REF] consider an explosion in the luggage compartment of an Aluminum-based fuselage, represented as a cantilever beam; their analysis neglects possible rigid-body motion of the whole structure. T. Kotzakolios and D. E. Vlachos [START_REF] Kotzakolios | Blast response of metal composite laminate fuselage structures using finite element modelling[END_REF] refer to Airbus A380 and introduce pressurization just as a static load on the skin.

In this paper, we propose a numerical procedure for evaluating the response of a fuselage subjected to an in-cabin explosion, with the aim of indicating a possible passive cabin protection. At variance of other approaches, our analysis includes gravity and pressurization loads at cruising altitude. We consider different volumes of air inside and outside cabin, different velocities of traveling shock waves, and changes in pressurization.

In Section 2 we describe fuselage's geometry and schemes for the pertinent design. Sections 3 and 4 deal with the representation of blast actions and the constitutive behaviour of Aluminum-based alloys, respectively. Section 5 describes possible passive protections based on Kevlar and polyurethane foams. We describe our numerical strategy in Section 6. Finally, in Section 7 we report simulations of in-cabin explosions and analyze the reliability of the proposed passive protection.

Fuselage geometry and design schemes

Fuselage's design refers commonly to three different schemes: truss, monocoque, and semimonocoque.

• The truss design, commonly belonging to the first generation of aircrafts, consists of steel tubes, welded together in a framework.

• The monocoque scheme refers to formers, frame assemblies, and bulkheads determining fuselage's shape.

• The semimonocoque is a modification of the latter design consisting of frame assemblies, bulkheads, and formers, supplemented by additional reinforcements, called longerons, which make the structure lightweighted and stiffer. Semimonocoque fuselages are usually made of aluminium alloys, although steel and titanium are used in high temperature regions. The fuselage is 4 m long and has a diameter of 3 m. Longerons and bulkheads appear in Figure 2, together with their sections. The floor consists of plates with 8 mm thickness, while the skin has 2 mm thickness. Tied contact pairs assure continuity between different parts. Al2024-T3 constitutes frames along the floor, longerons, and bulkheads.

Blast actions

Explosion produces a blast wave with high-pressure accompanying high-temperature expansion of gases. First, detonation induces a supersonic shock front.

With reference to a free-field explosion, Figure 3 shows a schematic representation of the hydrostatic overpressure P s = P -P o , i.e., the difference between the hydrostatic pressure P determined by the explosion and the ambient one, P o , as a function of the stand-off distance from the explosive.

The shock front is a discontinuity surface for the velocity field. Behind the wave front, a rarefaction wave propagates. Hydrostatic pressure and density decrease to values lesser than those in the ambient before the explosion. The hydrostatic overpressure, P s , at a point located at a distance R from the explosive decreases with both time t > t A (t A is the shock time arrival) and R. Generally, the time rate reduction is much greater than the spatial one. Figure 4 shows the schematic time variation of P s at a point. After a delay t A from detonation, the overpressure jumps suddenly from zero to P so . For t > t A , the overpressure decreases extremely fast until the instant t A + t o , the end of the so-called positive phase. Positive and negative impulses can be defined as the integrals of the hydrostatic overpressure along positive and negative phases, respectively:

i s = t A +to t A P s dt, i s-= t A +to+to- t A +to P s dt. (1)
Commonly, i s is greater than i s-. The duration of the negative phase is always longer than the positive one; also, |P so-| ≤ P so . The inequality i s >> i s-holds just in the near-field. The difference between the two impulses decays with distance.

The abrupt increase of the hydrostatic pressure within a blast wave can produce severe structural damage. When the primary shock meets a surface, it generates the so-called reflected overpressure, P r , which is characterized by a time-history similar to that of the hydrostatic overpressure. Commonly, the peak P ro of P r at time t A is much greater than the hydrostatic overpressure P so measured at the same point in the absence of any surface. The energy of a blast wave includes a contribution deriving from the hydrostatic overpressure, P s , and another one associated to from the dynamic pressure q := 1 2 ρu|u|, with ρ the density and u the velocity of gas particles. The hydrostatic overpressure is responsible for crushing of a surface. The dynamic pressure represents the amount of the kinetic energy of the stream-wise gas dynamics. It determines drag and possible lift forces. The ratio of this energy partition depends on the target geometry and its distance from the explosive charge: in the proximity of the source, P s and q attain almost equal values, while for increasing distances the dynamic pressure decreases until the blast wave decays to acoustic levels and the particle velocity becomes negligible. Furthermore, the drag loading induced blast wave depends on Mach and Reynolds numbers [START_REF] Dewey | Measurement of the physical properties of blast waves[END_REF][START_REF] Sochet | Blast effects. Physical properties of shock waves[END_REF][START_REF] Needham | Shock Wave and High Pressure Phenomena[END_REF]. We record different approaches to the simulation of a blast [START_REF] Remennikov | A review of methods for predicting bomb effects on buildings[END_REF][START_REF] Larcher | Simulation of the effects of an air blast wave[END_REF]. Here we refer to three phenomenological approaches: JWL, CONWEP, and TM5-1300. JWL stands for Jones, Wilkins, and Lee [START_REF] Jones | The detonation of solid explosives: the equilibrium conditions in the detonation wave-front and the adiabatic expansion of the products of detonation[END_REF][START_REF] Wilkins | The equation of state of PBX 9404 and LX 04-01[END_REF][START_REF] Lee | Adiabatic expansion of high explosives detonation products[END_REF]. It essentially rests on a state equation for the overpressure P s due to detonation:

P s = A 1 - ω ρ R 1 ρ 0 exp -R 1 ρ 0 ρ + B 1 - ω ρ R 2 ρ 0 exp -R 2 ρ 0 ρ + ω ρ e 0 . (2)
A, B, R 1 , R 2 , and ω are parameters depending upon the explosive, along with ρ 0 , its density, while ρ is the density of the detonation products and e 0 the internal energy per unit mass. All their values are selected to fit experimental results on a cylinder expansion test (see Table 2). Three different contributions determine P s : (i) the first exponential term accounts for high-pressure regime and high density of the detonation products; (ii) the second represents the intermediate pressure range; (iii) the last term is the constitutive equation of an ideal gas and defines the behaviour at low pressure or large expansion (see Fig. 5). The blast phenomenon can be modeled by taking advantage of JWL scheme through a multi-material approach, consisting of two domains: the explosive and the air, considered as an ideal gas with total pressure given by p

ρ ρ 0 A 1 ( R 1 ρ 0 ω ρ ( ( ρ ρ 0 ( R 1 exp B 1 ( R 2 ρ 0 ω ρ ( ( ρ ρ 0 ( R 2 exp ρ ω e0
+ p A = ρ RT, (3)
where p A and p are, respectively, ambient and current (fluid) pressures, R the product between Boltzmann's and Avogadro's constants, T the current absolute temperature. At the interface between air and explosive, shock waves occur and propagate [START_REF] Boyer | An experimental study of the explosion generated by a pressurized sphere[END_REF]. When we refer to JWL and a multi-material approach to the structural design, we need a fine discretization of both charge and fluid domain; and the latter can be wide.

Resource to empirical methods reduces computational costs. The two most commonly used empirical models rest upon related studies of the U.S. Army Corps of Engineers (USACE): the document [START_REF] Usace | Design and Analysis of Hardened Structures to Conventional Weapons Effects[END_REF] and the Technical Manual TM5-1300 [START_REF] Usace | Structures to Resist the Effects of Accidental Explosions[END_REF]. They contain the model CONWEP, completed by successive documents [START_REF] Usace | Structures to Resist the Effects of Accidental Explosions[END_REF]. The Joint Research Center of the European Union has produced in 2013 a technical report [START_REF] Karlos | Calculation of Blast Loads for Application to Structural Components[END_REF] substantially referring to these two last USACE documents and to another technicalrReport of the U.S. Army [START_REF] Kingery | Air blast parameters from TNT spherical air burst and hemispherical burst[END_REF].

Such empirical models do not allow a detailed reconstruction of the reflected waves, with the exception of ground effects (the case of hemispherical blast). However, depending upon geometry, the concentration of the reflected waves can give rise to local effects that can be greater than the original shock wave. In the case of vaulted structures at ground, limited laterally by walls, a localized shock wave produced by the reflected waves can hit the vault with an overpressure far greater than that produced directly by the original impinging shock wave [START_REF] Vannucci | A study on the simulation of blast actions on a monumental structure[END_REF], [START_REF] Vannucci | Cathédrales durables[END_REF], [START_REF] Masi | A study on the effects of an explosion in the Pantheon of Rome[END_REF].

In the case of aircrafts, due to the cylindrical symmetry, we can have prominent effects of focalisation. Moreover, the strong thermo-fluid-mechanics coupling requires an integrated methodology able to account for both stress in the structure and fast-dynamics explosion. For all these reasons, we adopt JWL for the numerical simulations.

Constitutive behaviour of some aluminum-based alloys 100

Among available proposals of constitutive structures for Aluminum-based alloys, we refer here to Johnson's and Cook's [START_REF] Johnson | A constitutive model and data for metals subjected to large strains, high strain rates and high temperatures[END_REF] scheme, which includes strain hardening, strain rate and thermal softening effects. We also add the elastic component of the constitutive behaviour.

Plastic behaviour

Johnson-Cook's scheme is a particular type of von Mises' plasticity. The yield stress, σ 0 , at zero plastic strain rate, is assumed to depend on the non-dimensional temperature, T , and the plastic component of the small strain, ε pl :

σ 0 = [JC 1 + JC 2 (ε pl ) n ] 1 -T m , (4)
where JC 1 is the yield strength, while JC 2 and n are, respectively, the hardening constant and an exponent to be determined on the basis of data fitting; m is the thermal softening exponent; also

T =    0 if T < T r (T -T r )/(T m -T r ) if T r ≤ T ≤ T m 1 if T > T m (5)
where T m and T r are melting and reference (ambient) temperatures, respectively. For T > T m , σ 0 = 0: the material melts. For nonzero plastic strain rate, εpl , we consider a correction to σ 0 , indicated by σ Y and defined by

σ Y = [JC 1 + JC 2 (ε pl ) n ] 1 + JC 3 ln εpl ε0 1 -T m , ( 6
)
where ε0 is a reference strain rate, adopted for normalization purposes.

Damage

We associate damage to plastic strain. We define its increment, δ, by

δ = k=K k=0 ∆ε pl,k ε pl,f , (7)
where ∆ε pl is an increment of the equivalent plastic strain, ε pl,f the plastic strain at failure, and K the total number of steps in the analysis. By assumption, the strain at failure is given by

ε pl,f = d 1 + d 2 exp d 3 p q 1 + d 4 ln εpl ε0 1 + d 5 T , (8)
where p indicates pressure, q von Mises yield stress, and d i , i = 1, . . . , 5, are failure parameters. Tables 4 and5 show the pertinent parameters for Al7075-T6 and Al2024-T3.

Passive protections 110

For their high toughness, we often make use of aramidic fibers in lightweight composites adopted in aerospace engineering. Kevlar failure emerges from a sequence of ruptures of sub-fibers. The resulting behaviour is ductile, with pertinent dissipation of energy in time. For this reason, we consider here a protection based on Kevlar fabric. We imagine a cabin load device (CLD) constructed by inserting Kevlar in the cabin, tied with longerons and bulkheads, upholstering the interiors. During an explosion, while experiencing large strain, CLD may protect the skin and prevent from fuselage perforation, as we shall see in Section 7. We can increase the protection by using foam in cavities of the fuselage such as the spaces between longerons and the interior of the floor. Foams, as energy absorbing materials, allow to mitigate the stresses due to impact, being characterised, at the same time, by low densities. Polyurethane foams are excellent candidates for blast protections.

Kevlar

We assume for Kevlar fabric an orthotropic linear elastic behaviour up to failure. Also, we consider a plane stress state (σ 13 = σ 23 = σ 33 = 0). With ε the total elastic strain, and C the fourth-rank elasticity tensor, for the Cauchy stress tensor, σ, we write

σ = C ε (9)
With E i Young moduli, ν ij Poisson ratios, and G ij shear moduli, i, j = 1, 2, 3, in Voigt's notation, with

γ ij = 2 ij , relation (9) writes    ε 11 ε 22 γ 12    =   1/E 1 -ν 21 /E 2 0 -ν 12 /E 1 1/E 2 0 0 0 1/G 12      σ 11 σ 22 σ 12    (10)
and

ε 33 = - ν 13 E 1 σ 11 - ν 23 E 2 σ 22 = 0 (11)
The elastic parameters, collected in Table 6, refer to a woven fabric Kevlar-29 (style 735), with a polyvinylbutyral-phenolic matrix, with 18% fibers. We refer to laminated plates fabricated at a room temperature curing polyester matrix by wet lay-up techniques, with a square-weave fabric, i.e., equal number of fibers along the two orthogonal directions woven in a 2-harness satin style, where each fiber bundle loops alternately over and under two bundles in the orthogonal direction. The thickness of the fabric is approximately 0.61 mm, with an areal density of 0.491 kg m -2 . The resin matrix is a low viscosity thermosetting polyester resin (1.22 mm total thickness) [START_REF] Zhu | Penetration of laminated Kevlar by projectiles -I. Experimental investigation[END_REF].

We consider Hashin's failure criterion for fiber-reinforced materials [START_REF] Hashin | Failure Criteria for Unidirectional Fiber Composites[END_REF], so we account for four different damage mechanisms: fiber tension and compression, matrix tension and compression. The initiation criterion relies on a set of damage variables, each one corresponding to a damage mechanism:

d t f = σ 11 X t 2 + σ 12 S l 2 , ( 12
)
d c f = σ 11 X c 2 , ( 13
)
d t m = σ 22 Y t 2 + σ 12 S l 2 , (14)
d c m = σ 22 2S t 2 + Y c 2S t 2 -1 σ 22 Y c + σ 12 S l 2 , (15)
where X t and X c are, respectively, longitudinal tensile and compressive strengths, Y t and Y c transverse tensile and compressive strengths, S l and S t in-plane and transversal shear strengths, σ 11 , σ 22 , and σ 12 the components of the effective stress tensor defined by

σ = Mσ, ( 16
)
where M is the damage operator, which we assume to be

M =   1/(1 -d f ) 0 0 0 1/(1 -d m ) 0 0 0 1/(1 -d s )   . (17)
d f , d m and d s are damage variables associated with fiber, matrix and shear modes, respectively, namely

d f = d t f if σ 11 ≥ 0 d c f if σ 11 < 0 ( 18
)
d m = d t m if σ 22 ≥ 0 d c m if σ 22 < 0 ( 19
)
d s = 1 -1 -d t f 1 -d c f 1 -d t m (1 -d c m ) (20)
Table 7 collects the parameters of Hashin's failure criterion, together with the fracture energy G f [START_REF] Wang | Mechanical properties of fiber glass and Kevlar woven fabric reinforced composites[END_REF], assumed to be equally distributed in longitudinal, transverse and shear directions.

Polyurethane foam

As usual, we consider polyurethane as a visco-elastic material with null Poisson ratio and a density of 240 kg/m 3 . Data about the constitutive coefficients derive from uniaxial tension and compression tests, at different strain rates [START_REF] Goods | Mechanical properties and energy absorption characteristics of a polyurethane foam[END_REF]. We assume a maximal value of 3.8 MPa for the maximum principal tensile stress. Furthermore, we consider strain rate-dependent behaviour under compression.

Numerical procedure

We simulate an in-cabin explosion by using a CEL approach: the fuselage is immersed into air; we use JWL to model the blast.

In CEL simulations, ABAQUS/Explicit takes into account the Eulerian fluid domain through the socalled volume-of-fluid-method: an additional variable within each element, the Eulerian volume fraction (EVF), describes the flow through the mesh. In our case, this approach allows us to compute not only the propagation of shock waves, but also the diffusion of the explosive inside the air domain. We describe fluid-structure interaction by using a general contact algorithm, with a null interface friction coefficient and a penalty method.

We focus our attention just on a portion of the fuselage. For the sake of simplicity, we do not include seats, passengers, and baggage, so that we do not consider the projectile-like behaviour of fragments of bodies invested by an explosion. We do not account for bolts, rivets, and points of adhesive bonding, which are possible source of stress concentration, influencing the damage pattern. This is a limitation of our numerical analyses. However, our aim is to present a procedure, showing its potentialities, rather than designing a specific aircraft.

For the structure under analysis, we account for the stress state at cruising altitude, before explosion. Our methodology consists of two steps:

• First, we develop a quasi-static analysis by applying increments of gravity and a distributed load equivalent to lift, on the lower half of the fuselage. Figure 8 shows boundary conditions, which allow a rigid body motion as soon as vertical loads become unbalanced.

• The results are the initial state of a subsequent simulation in which we include pressurization by means of a pressure gradient between internal and external volumes of air. We consider the aircraft to be at 10000 m altitude, while the equivalent effective cabin is at 2000 m.

Figure 9 is a pertinent scheme. We use the stress state obtained from the second analysis as initial state of the simulations of an internal explosion at cruising altitude.

To account for blast heat release, we adopt a Coupled-Temperature-Displacement (CTD) analysis. The equilibrium state achieved at the end of the second step is the initial state of a CEL-CTD explicit analysis of blast waves on the structure. Table 2 collects constitutive parameters for explosive and air.

The choice of the mesh fineness derives from a convergence analysis. We have considered six different meshes for the air domain: from the coarsest, M 1 , with 0.58 × 10 6 DOF, to the finest, M 6 , with 10.54 × 10 6 DOF. For each mesh, we have computed the blast overpressure due to an explosive charge placed in the middle of the fuselage length, identified by plane 1. Figure 10 indicates charge and gauge points (further details are in Figure 11).

Figure 12 shows the relative error evaluated for meshes M 1 , . . . , M 5 with respect to mesh M 6 , and defined by

err = 100 p Mi -p M6 p M6 , (21)
where p Mi is the peak overpressure in the gauge point evaluated for grid M i , while p M6 the corresponding value in the reference mesh M 6 . In showing results, we refer to the mesh M 5 , with 5.72 × 10 6 DOF. Figure 12 shows also a comparison between the overpressure obtained by using JWL and the analogous result emerging from an empirical model based on Kingery's and Bulmash's relations, complemented by Frielander's equation [START_REF] Sachs | The dependence of blast on ambient pressure and temperature[END_REF]. In such a scheme, Z Rp 1/3 0 ρe0 plays a role. R is the geometric distance from the explosive, p 0 the ambient pressure, ρ the density, e 0 already introduced in relation 2 [START_REF] Sachs | The dependence of blast on ambient pressure and temperature[END_REF]. Moreover, we have considered the impulse reduction due to the angle of incidence of the wave over the structure. By comparing the two analyses, we record for JWL a lower pressure peak, a clearing-type effect, and an almost negligible negative phase. We do not consider confined explosions, because the fuselage may be perforated, reducing the confinement (see also the remarks in Appendix).

Mesh sensitivity of the structural domain

For the fuselage we have considered four-nodes thermally coupled shell elements, supported by the explicit solver of ABAQUS. A numerical convergence analysis corroborates the choice. As already pointed out, first we consider gravity and lift loads. Then, the equilibrium state enters a second analysis in which we consider the fluid domain, discretized with the selected mesh M 5 , and we apply the pressurization load.

The convergence analysis rests on explicit schemes for nine different structural meshes MS i , i = 1, . . . , 9, with characteristics in Table 8.

We have evaluated convergence as for the air domain, computing for each mesh MS i , with i = 1, . . . , 8, the errors ∆ε pertaining to the evaluation of the strain ε, relatively to a reference ε r , calculated for the finest mesh, MS 9 . We define

∆ε = ε -ε r ε r , (22)
for any mesh MS i . The strain refers to the point indicated in Figure 13. The same figure shows the result. Consequently, we have selected the mesh MS 8 for the final simulations. A detail of the structural mesh MS 8 appears in Figure 14.

The explosive has an initial square shape and has been discretized by using elements of 3 mm. We use parallel computing for decreasing significantly computational time. All the simulations exploited a 24-core workstation.

mid-length plane 330 300 gauge point mid-length section

Numerical results

First, we have simulated a blast, produced by an in-cabin explosion of 0.5 kg TNT equivalent, by using a CEL-CTD analysis. Then, we have designed a CLD (Fig. 6) by suing the results. The stable time step of the explicit schemes, mostly afflicted by the small elements size (see [START_REF] Abaqus | ABAQUS Documentation[END_REF]), varies from 2.9 × 10 -7 s to 5.7 × 10 -11 s. Therefore, the computational time required for running a short time period (3.5 milliseconds) analysis is about 35 hours. The blast is really fast and so is fuselage's reaction, which can be evaluated already in the first 3 milliseconds. The final structural damage is in Figure 16, while the progressive response of bulkheads, longerons and skin to the blast is in Figure 17. We do not consider cutouts along which longerons can run. The results confirm that a small amount of charge, capable just of piercing the skin at sufficiently high altitude, can be the origin of a disaster. The principal shock front produces a hole in the skin enlarged by the subsequent decompression until global failure. The overpressure destroying the skin around the charge may break also the floor (see Fig. 15). In this way it determines a secondary shock, which propagates inside the cargo hold and opens another hole in the bottom half of the skin.

In-cabin explosion simulations

Figures 18,19, and 20 refer to the analogous analysis developed foreseeing the sole Klevar-based protection.

If we use just Kevlar, in protecting the fuselage, we have an increment in weight of about 1.2 per cent, while the presence of foam increase weight up to 38 per cent.

Assessing the reliability of the protective CLD

Figure 22 shows the evolution of the pressure field in the air domain considering the protective CLD: we get mitigation of blast waves, thanks to both materials and the way the protections are inserted. By comparing Figures 22 and15, we can see the advantages of using a Kevlar fabric as a result of the smooth evolution of the principal shock front and reflected waves with a lower magnitude. Indeed, such a fabric behaves like a second skin, much less stiff and tougher than the one made of Al7075-T6. Its main effect is a reduction of the magnitude of reflected blast waves.

Figures 23 and24 show the structural damage and the progressive response of fuselage and CLD, respectively. Skin can withstand an internal explosion of 0.5 kg, being not punctured, and floor is only partially damaged, avoiding the failure of the bottom half skin, due to venting, as shown in Figure 22, where any blast wave penetrates inside the cargo hold.

Figure 25 displays the time-evolution of the energy dissipated along the damage process in the three cases here studied: a) without protections, b) with Kevlar, and c) with the protective CLD. A Kevlar-based protection (case b) results in a reduction of the amount of structural damage equal to ∼ 79% with respect to case a. The protective CLD reduces further the zone of the fuselage affected by damage by ∼ 46% of the one computed without protections.

We remark that in performing our analyses by assuming just static pressure load, we do not record damage. The circumstance underlines the difference between our approach and those based on static pressure load.

Further remarks

Numerical simulations dealing with the behaviour of fuselages under blasts are usually based on the application of pressurization load as a distributed action on the skin independently of the pressures of the surrounding volumes of air. We investigate this assumption in Appendix.

The simulations presented here allow us to capture the fast-dynamics due to an internal explosion. The difference between numerical results and reality requires reliable experimental tests. However, our results furnish indication addressing towards a reasonably profitable design.

We have used just a portion of a fuselage free of cabin furnitures, passengers, and baggage. Also, in our analyses we did not include after-burning, conduction, and chemical effects. These choices and those already mentioned are limitations of our work. However, our aim is just to indicate a nonstandard methodology of analysis, showing its potentialities, rather than designing a specific airplane.

We record (ground) tests on Kevlar-reinforced cabins as a blast mitigation measure [START_REF] Rediess | Aircraft blast mitigation[END_REF], [START_REF]Aircraft Hardening Research Programme Final Overview Report[END_REF]. Our analyses at a cruising altitude corroborate that choice. However, the methodology that we propose can be used to analyze other possible mitigation measures. [START_REF] Corona | An Evaluation of the Johnson-Cook Model to Simulate Puncture of 7075 Aluminum Plates[END_REF] and [START_REF] Kray | Failure modeling of titanium 6Al-4V and aluminum 2024-T3 woth the Johnson-Cook material model[END_REF]). Table 6: Material parameters for Kevlar fabric (according to [START_REF] Zhu | Penetration of laminated Kevlar by projectiles -I. Experimental investigation[END_REF] and [START_REF] Cheeseman | Ballistic impact into fabric and compliant composite laminates[END_REF]). (see Fig. 29). In any case, however, the numerical results show that the procedure of defining two distinct volumes of air can take into account the abrupt changes in the pressure field (thus in the pressurization load). At variance, we would neglect a physically significant physical contribution (abrupt pressure variation).

Material

100×100 mm hole

Figure 1 :

 1 Figure 1: Truss design (left) and semimonocoque design (right).

Figure 2 :

 2 Figure 2: Fuselage, longerons and frame assemblies and related sections.

Figure 3 :

 3 Figure 3: Spatial distribution of the overpressure.

Figure 4 :

 4 Figure 4: Scheme of the time variation of the overpressure due to a blast.

Figure 5 :

 5 Figure 5: Overpressure inside an explosive material, 2,4,6-trinitrotoluene, at varying of the density of the detonation products.

Figure 6 :

 6 Figure 6: CLD and covering in polyurethane foam.

  130

Figure 7 :

 7 Figure 7: Behaviour of polyurethane foam in compression and tension from [41].

Figure 10 :Figure 11 :

 1011 Figure 10: The air domain: overall view (left), and location of the gauge and source points (right).

Figure 12 :

 12 Figure 12: Relative error, eq. (21), for the gauge point as function of the number of DOF (left). Comparison of the overpressure time history at the gauge point between the numerical simulation with mesh M 5 and empirical model based on the relationships of Kingery-Bulmash (right).

Figure 13 :

 13 Figure 13: Relative error ∆ε versus DOF, with respect to a point of the skin.

Figure 14 :

 14 Figure 14: Detail of the structural mesh MS 8 .

Figures 15 ,

 15 Figures 15, 16 and 17 report the results. The evolution of the pressure field thorugh the air domain is in Figure 15.

Figure 15 :

 15 Figure 15: The evolution of the pressure field through the air domain.

Figure 16 :

 16 Figure 16: Structural damage at the end of the calculations.

Figure 19 :

 19 Figure 19: The evolution of the pressure field through the air domain, with the liner in Kevlar.

Figure 20 :

 20 Figure 20: Structural damage at the end of the calculations, with the liner in Kevlar.

Figure 21 :

 21 Figure 21: Simulation of an internal blast with the liner in Kevlar: external (left), internal (center ), and front views (right).

Figure 22 :

 22 Figure 22: The evolution of the pressure field through the air domain, with CLD.

Figure 23 :

 23 Figure 23: Structural damage at the end of the calculations, with CLD.

Figure 24 :

 24 Figure 24: Simulation of an internal blast with CLD: external (left), internal (center ), and front views (right).

Figure 25 :

 25 Figure 25: Evolution of the (normalized) energy dissipated through damage in the three cases: without protections, with Kevlar, and with the protective CLD.

d 1 d 2 d 3 d 4 d 5

 5 Elongation at break[%] Al7075-T6 0.005 0.34 -1.5 -0.039 8 22 Al2024-T3 0.13 0.13 -1.5 0.011 0 18

Figure 26 :

 26 Figure 26: Location of the square hole with zero strength.

Figure 28 :

 28 Figure 28: Structural response after the removal of a 100 × 100 square hole from the skin.

  At the instant t A + t o , the so-called negative phase starts. It pertains to the rarefaction wave, triggered by the expansion of the detonation products: P s decreases to negative values and asymptotically approaches zero after t A + t o + t o-.

	70
		P so
	s Overpressure P	i s
	t A	t A +t o	i s-	t A +t o +t o-
	0			time
		P so-
		positive	negative
		phase	phase

Table 1 :

 1 Non-accidental explosions in airplanes, from[9].

	Date Flight	Description	Casualities
	Boeing 247D, United Air Lines	bomb made of nytroglicerin	7
		placed in the baggage compartment
	Douglas C-47-DL,	bomb made of dynamite	23
	Canadian Pacific Airlines	placed in the baggage compartment
	Douglas DC-6B, United Air Lines	bomb made of dynamite	44
	Boeing 707-124, Continental Air Lines explosive device inside passenger cabin	45
	Douglas RD4-1, Aden Airlines	explosive device inside passenger cabin	30
	DH-106 Comet 4,	high explosive device	66
	British European Airways	within the cabin under seats
	Convair CV-990-30A-6 Conorado,	bomb in the baggage compartment	47
	Swissair
	Douglas DC-8-43, Cubana de Aviación explosive device at the rear of the cabin	73
	Boeing 720-023B, Middle East Airlines bomb in the baggage compartment	81
	Boeing 747-121, Pan Am	bomb placed under a seat cushion	1
	Boeing 747-237B, Air India	high explosive device inside	329
		the cargo compartment
	Boeing 727-231, Trans World Airlines	explosive device in the cabin	4
	Boeing 737-270C, Iraqi Airways	two hand grenades in the cabin	63
	Boeing 707-3B5C, Korean Air	liquid explosives concealed	115
		as liquor bottles
	Boeing 747-121, Pan Am	high-explosive device in the cabin	270
	McDonnel Douglas DC-10-30,	high explosive device	170
	Union de Transport Aériens	in the cargo hold
	Boeing 727-21, Avianca Airlines	explosive near the fuel tank	110
	Boeing 767, American Airlines	plastic explosive concealed	0
		within shoes
	Tupolev Tu-134-3, Volga-Avia Express high explosive	90
	Tupolev Tu-154B-2, Siberia Airlines	devices
	Airbus A32-231, Metrojet	1 kg of TNT	224
	Airbus A321-111, Daallo Airlines	explosive device concealed within	1
		a laptop computer

Table 2 :

 2 Physical parameters for the state equations for air an explosive, after[START_REF] Rogers | Thermodynamic and Transport Properties of Fluids[END_REF] and[START_REF] Lee | JWL equation of state coefficient for high explosives[END_REF].

	AIR	ρ	Temperature Gas constant Specific heat
		[kg/m 3 ]	[K]	[J/(kg K)]	[W/(m K)]
		1.225	288.2	287.058	717.6
	TNT	ρ 0	v D	A	B	ω	R 1	R 2	e 0
		[kg/m 3 ]	[m/s]	[MPa]	[MPa]				[kJ/kg]
		1630	6930	373770	3747.1	0.35 4.15 0.9	3680

Table 3 :

 3 Material elastic and thermal parameters.

	Material	ρ	E	ν	Specific	Thermal
					heat	conductivity
		[kg/m 3 ] [GPa]		[J/(kg K)]	[W/(m K)]
	Al7075-T6	2810	71.7	0.33	848	130
	Al2024-T3	2770	73.1	0.33	863	121

Table 4 :

 4 Material parameters for Johnson-Cook's plasticity (according to[START_REF] Corona | An Evaluation of the Johnson-Cook Model to Simulate Puncture of 7075 Aluminum Plates[END_REF] and[START_REF] Kray | Failure modeling of titanium 6Al-4V and aluminum 2024-T3 woth the Johnson-Cook material model[END_REF]).

	Material	JC 1	JC 2	n	m	JC 3	T r	T m	ε0
		[MPa] [MPa]				[K]	[K]	[1/s]
	Al7075-T6	517	405	0.41 1.1 0.0075 292.2 750 0.000161
	Al2024-T3	265	426	0.34 1.0 0.015 293.2 775	1.0

Table 5 :

 5 Material parameters for Johnson-Cook's damage (according to

Table 7 :

 7 Material parameters for Hashin's failure criterion[START_REF] Cheeseman | Ballistic impact into fabric and compliant composite laminates[END_REF].

	X t	X c	Y t	Y c	S t	S l	G f
	[MPa] [MPa] [MPa] [MPa] [MPa] [MPa] [kJ/m 2 ]
	555	555	34.5	34.5	898	77	1.64

Table 8 :

 8 Characteristics of the studied structural meshes.

	Mesh Average element size [mm] Number of elements Number of DOF [mln]
	MS 1	10	9191	7.23 × 10 -2
	MS 2	8	12520	1.01 × 10 -1
	MS 3	6	21329	1.73 × 10 -1
	MS 4	5	27971	2.24 × 10 -1
	MS 5	4	40813	3.22 × 10 -1
	MS 6	3	74611	5.87 × 10 -1
	MS 7	2	159764	1.26
	MS 8	1.5	266811	2.10
	MS 9	1	597159	4.7

Mesh sensitivity of the fluid domain

The air domain is shown in Figures 10 and11. It is decomposed into 1.156 × 10 6 8-nodes hexaedral volume elements, for 5.78 × 10 6 degrees of freedom (DOF).

Boundary conditions consist on non-reflecting outflow, applied to the whole boundary. We subdivide the air domain into two different volumes: inside and outside the fuselage. In the latter case, the domain corresponds to a cylinder with a diameter of 4 m.

Appendix

We analyze here the numerical procedure used in the above numerical simulations and presented in Section 6. We pay attention to the use of two distinct volumes of air to represent the pressurization load 260 in terms of a pressure gradient. We aim at a reasonable characterization of the pressure evolution on the cabin, which may be decompressed, due to the explosion-induced damage. Hence, we simulate a sudden loss of internal pressure.

We proceed as follows: First, we compute the stress distribution in the fuselage, by applying gravity, lift, and pressurization loads, along the path described in Section 6. Then, we impose a 100 × 100 mm square hole over the skin by equating to zero the strength of the corresponding finite elements (see Fig. 26).

The air flows outside the hole, due to the pressure gradient, at high velocity. The material around the hole undergoes plastic strain; damage increases. Figure 27 shows the evolution of the pressure field immediately after the puncturing, while the structural response is displayed in Figure 28.

The simulations do not describe the whole decompression event; calculations end, in fact, after 30 ms