Vitor Paisante
email: paisante@dcc.ufmg.br

Maroua Maalej
email: maroua.maalej@ens-lyon.fr

Leonardo Barbosa

Laure Gonnord
email: laure.gonnord@ens-lyon.fr

Fernando Magno
email: fernando@dcc.ufmg.br

Quintao Pereira

Symbolic Range Analysis of Pointers *

Keywords: Alias analysis, range analysis, speed, precision

Alias analysis is one of the most fundamental techniques that compilers use to optimize languages with pointers. However, in spite of all the attention that this topic has received, the current state-of-the-art approaches inside compilers still face challenges regarding precision and speed. In particular, pointer arithmetic, a key feature in C and C++, is yet to be handled satisfactorily. This paper presents a new alias analysis algorithm to solve this problem. The key insight of our approach is to combine alias analysis with symbolic range analysis. This combination lets us disambiguate fields within arrays and structs, effectively achieving more precision than traditional algorithms. To validate our technique, we have implemented it on top of the LLVM compiler. Tests on a vast suite of benchmarks show that we can disambiguate several kinds of C idioms that current state-of-the-art analyses cannot deal with. In particular, we can disambiguate 1.35x more queries than the alias analysis currently available in LLVM. Furthermore, our analysis is very fast: we can go over one million assembly instructions in 10 seconds.

Introduction

Pointer analysis is one of the most fundamental compiler technologies. This analysis lets the compiler distinguish one memory location from others; hence, it provides the necessary information to transform code that manipulates memory. Given this importance, it comes as no surprise that pointer analysis has been one of the most researched topics within the field of compiler construction [START_REF] Hind | Pointer analysis: Haven't we solved this problem yet?[END_REF]. This research has contributed to make the present algorithms more precise [START_REF] Hardekopf | The ant and the grasshopper: fast and accurate pointer analysis for millions of lines of code[END_REF][START_REF] Zhang | Efficient subcubic alias analysis for C[END_REF], and faster [START_REF] Hardekopf | Flow-sensitive pointer analysis for millions of lines of code[END_REF][START_REF] Shang | On-demand dynamic symmarybased points-to analysis[END_REF]. Nevertheless, one particular feature of imperative programming languages remains to be handled satisfactorily by the current state-of-the-art approaches: the disambiguation of pointer intervals.

Mainstream compilers still struggle to distinguish intervals within the same array. In other words, state-of-the-art pointer analyses often fail to disambiguate regions addressed from a common base pointer via different offsets, as explained by Yong and Horwitz [START_REF] Yong | Pointer-range analysis[END_REF]. Field-sensitive pointer analysis, provide a partial solution to this problem. These analyses can distinguish different fields within a record, such as a struct in C [START_REF] Pearce | Efficient fieldsensitive pointer analysis for C[END_REF], or a class in Java [START_REF] Yan | Demand-driven contextsensitive alias analysis for java[END_REF]. However, they rely on syntax that is usually absent in the low level program representations adopted by compilers. Shape analyses [START_REF] Jones | A flexible approach to interprocedural data flow analysis and programs with recursive data structures[END_REF][START_REF] Sagiv | Solving shape-analysis problems in languages with destructive updating[END_REF] can disambiguate subparts of data-structures such as arrays, yet their scalability remains an issue to be solved. Consequently, many compiler optimizations, such as loop transformations, tiling, fission, skewing and interchanging [START_REF] Wolfe | High Performance Compilers for Parallel Computing[END_REF]Ch.09], are very limited in practice. Therefore, we claim that, to reach their full potential, compilers need to be provided with more effective alias analyses.

This paper describes such an analysis. We introduce an abstract domain that associates pointers with symbolic ranges. In other words, for each pointer p we conservatively estimate the range of memory slots that can be addressed as an offset of p. We let GR(p) be the global abstract address set set associated with pointer p, such that if loc i + [l, u] ∈ GR(p), then p may dereference any address from @(loc i ) + l to @(loc i ) + u, where loc i is a program site that contains a memory allocation call, and @(loc i ) is the actual return address of the malloc at runtime. We let {l, u} be two symbols defined within the program code. Like the vast majority of pointer analyses available in the compiler literature, from Andersen's work [START_REF] Andersen | Program Analysis and Specialization for the C Programming Language[END_REF] to the more recent technique of Zhang et al. [START_REF] Zhang | Efficient subcubic alias analysis for C[END_REF], our method is correct if the underlying program is also correct. In other words, our results are sound with respect to the semantics of the program if this program has no undefined behavior, such as out-of-bounds accesses.

The key insight of this paper is the combination of pointer analysis with range analysis on the symbolic interval lattice. In a symbolic range analysis, ranges are defined as expressions of the program symbols, a symbol being either a constant or the name of a variable. There exist many approaches to symbolic range analyses in the literature [START_REF] Blume | Symbolic range propagation[END_REF][START_REF] Nazaré | Validation of memory accesses through symbolic analyses[END_REF][START_REF] Rugina | Symbolic bounds analysis of pointers, array indices, and accessed memory regions[END_REF]. The algorithms that we present in this paper do not depend on any particular implementation. Nevertheless, the more precise the range analysis that we use, the more precise the analysis facts that we produce. In this work we have adopted the symbolic range analysis proposed in 1994 by William Blume and Rudolf Eigenmann [START_REF] Blume | Symbolic range propagation[END_REF].

To validate our ideas, we have implemented them in the LLVM compilation infra-structure [START_REF] Lattner | LLVM: A compilation framework for lifelong program analysis & transformation[END_REF]. We have tested our pointer analysis onto three different benchmarks used in previous work related to pointer disambiguation: Prolangs [START_REF] Ryder | A schema for interprocedural modification side-effect analysis with pointer aliasing[END_REF], PtrDist [START_REF] Zhao | Dynamic memory optimization using pool allocation and prefetching[END_REF] and MallocBench [START_REF] Grunwald | Improving the cache locality of memory allocation[END_REF]. As we show in Section 4, our analysis is linear on the size of programs. It can go over one-million assembly instructions in approximately 10 seconds. Furthermore, we can disambiguate 1.35x more queries than the alias analysis currently available in LLVM.

Overview

We have two different ways to answer the following question: "do pointers tmp i and tmp j alias?" These tests are called global and local. In this section, we will use two different examples to illustrate situations in which each query is more effective. These distinct strategies are complementary: one is not a superset of the other. particular example, messages have two parts: an identifier, which is stored in the beginning of the array, and a payload, which is stored right after. The loops in lines 5-8 and 9-12 fill up each of these parts with data. If a compiler can prove that the stores at lines 6 and 10 are always independent, then it can perform optimizations that would not be possible otherwise. For instance, it can parallelize the loops, or switch them, or merge them into a single body.

Global pointer disambiguation.

No alias analysis currently available in either gcc or LLVM is able to disambiguate the stores at lines 6 and 10. These analyses are limited because they do not contain range information. The range interval [l, u] associated with a variable i is an estimate of the lowest (l) and highest (u) values that i can assume throughout the execution of the program. In this paper, we propose an alias analysis that solves this problem. To achieve this goal, we couple this alias analysis with range analysis on symbolic intervals [START_REF] Blume | Symbolic range propagation[END_REF]. Thus, we will say that the store at line 6 might modify any address from p + 0 to p + N -1, and that the store at line 10 might write on any address from p + N to p + N + strlen(m) -1. For this purpose, we will use an abstract address that encodes the actual value(s) of p inside the prepare function. These memory addresses are depicted in Figure 2, where each represents a memory slot.

Whole program analysis reveals that there are two candidate locations that any pointer in the program may refer to. These locations have been created at lines 17 and 18 of Figure 1, and we represent them abstractly as loc 17 and loc 18 . These names are unique across the entire program. After running our analysis, we find out that the abstract state (GR) of i at line 6 is GR(i n.6 ) = {loc 17 + [0, N -1], and that the abstract state of i at line 10 is GR(i n.10 ) = {loc 17 + [N, N + strlen(m) -1]}. Given that these two abstract ranges do not intersect, we know that the two stores update always different locations. We call this check the global disambiguation criterion.

Local pointer disambiguation. Figure 3 shows a program in which the simple intersection of ranges would not let us disambiguate pointers tmp 0 and tmp 1 . After solving global range analysis for that program, we have that GR(tmp 0 ) = {loc 0 + [0, N + 1]} and that GR(tmp 1 ) = {loc 0 + [1, N + 2]}, where loc 0 defines the abstract address of the function parameter p. The intersection of these ranges is non-empty for N ≥ 1. Thus, the global check that we have used to disambiguate locations in Figure 1 does not work in Figure 3. Notwithstanding this fact, we know that tmp 0 and tmp 1 will never point to a common location. In fact, these pointers constitute different offsets from the same base address. To deal with this imprecision of the global check, we will be also discussing a local disambiguation criterion. In this case, we rename every pointer p that is alive at the beginning of a single entry region to a fresh name new p . Whereas we use the global test for pointers in different regions, the local test is applied onto pointers within the same single entry region. After renaming, we update the table of pointer pairs, so that LR(new p ) = loc new + [0, 0], regardless of the old ranges assigned to the original pointer p. In Figure 4 we would have that LR(tmp 2 ) = loc new + [0, 0] and LR(tmp 3 ) = loc new + [START_REF] Aho | Compilers: Principles, Techniques, and Tools[END_REF][START_REF] Aho | Compilers: Principles, Techniques, and Tools[END_REF], where tmp 2 is the name of the address new p [0], and tmp 3 is the name of the address new p [START_REF] Aho | Compilers: Principles, Techniques, and Tools[END_REF]. This new binding of intervals to pointers gives us empty intersections between similar locations in LR(tmp 2 ) and LR(tmp 3 ). Consequently, the local check is able to distinguish addresses referenced by tmp 2 and tmp 3 .

Combining Range and Pointer Analyses

We perform our pointer analysis in several steps. Figure 5 shows how each of these phases relates to each other. Our final product is a function that, given two pointers, p 0 and p 1 , tells if they may point to overlapping areas or not. An invocation of this function is called a query. We use an offthe-shelf symbolic range analysis, e.g., à la Blume [START_REF] Blume | Symbolic range propagation[END_REF], to bootstrap our pointer analysis. By inferring the symbolic ranges of pointers, we have two alias tests: the global and the local approach. In the rest of this section we describe each one of these contributions.

A Core Language

We solve range analysis through abstract interpretation. To explain how we abstract each instruction in our intermediate representation, we shall use the language seen in Figure 6; henceforth, we shall call this syntax our core language. We shall be working on programs in Extended Static Single Assignment (e-SSA) form [START_REF] Bodik | ABCD: eliminating array bounds checks on demand[END_REF]. E-SSA form is a flavor of Static Single Assignment (SSA) [START_REF] Cytron | Efficiently computing static single assignment form and the control dependence graph[END_REF] form, with variable renaming after inequalities. Thus, our core language contains φfunctions to ensure the single definition (SSA) property, and intersections to rename variables after conditionals. We assume that φ-functions have only two arguments. Generalizing this notation to n-ary functions is immediate.

= malloc(i 0 ) -Free memory | p 0 = free(p 1 ) -Pointer plus int | p 0 = p 1 + i 0 -Pointer plus const | p 0 = p 1 + c 0 -Bound intersection | p 0 = p 1 ∩ [l, u] -Load into pointer | p 0 = * p 1 -Store from pointer | * p 0 = p 1 -φ-function | p 0 = φ(p 1 : 1 , p 2 : 2 ) -Branch if not zero | bnz(v, ) -Unconditional jump | jump( ) Figure 6
. The syntax of our language of pointers. Figure 7 shows the control flow graph of the program seen in Figure 1. The implementation of the analysis that we shall present in this paper is interprocedural, albeit not contextsensitive. To achieve interprocedurality, we associate actual parameters with formal parameters of functions. In the example of Figure 1, pointer b -an actual parameter -is linked with p -a formal parameter -through a φ-function.

Z = atoi(...) b = malloc(Z) "loc 0 " s = malloc(...) "loc 1 " p = b N = Z m 0 = s i 0 = p e = p + N i 1 = ϕ(i 0 , i 3 ) (i 1 < e)? i 2 = i 1 ∩ [-∞, e -1] *i 2 = 0 t 0 = i 2 + 1 *t 0 = 0xFF i 3 = i 2 + 1 i 4 = i 1 ∩ [e, +∞] f = e + strlen(m 0 ) i 5 = ϕ(i 4 , i 7 ) m 1 = ϕ(m 0 , m 2 ) (i 5 < f)? i 6 = i 5 ∩ [-∞, f -1] *i 6 = *m 1 m 2 = m 1 + 1 i 7 = i 6 + 1
The e-SSA format lets us implement our analysis sparsely, e.g,. we can assign information directly to variables, instead of pairs of variables and program points. As demonstrated by Choi et al. [START_REF] Choi | Automatic construction of sparse data flow evaluation graphs[END_REF], the main advantage of a sparse analysis is efficiency: the product of the analysis -the information that is bound to each variable -requires O(N ) space, where N is the number of variable names in the program. Furthermore, as we shall explain in the rest of this section, our analysis can be computed in O(N ) time.

Key to these good properties is the fact that we create new variable names at each program point where our analysis can infer new information. This knowledge appears due to memory allocation (malloc), deallocation (free), pointer arithmetic, intersections and φ-functions. Each of these instructions defines new variables, whose names are associated with information. For instance, the instruction p 0 = free(p 1 ) copies p 1 to p 0 , and binds p 0 to a memory chunk of size 0. As we will show in Section 3.4, our abstract interpreter associates with p 0 a new abstract state which indicates that p 0 is not a valid reference to any location.

Program Locations.

Our analysis binds variable names to sets of locations and ranges. We denote the set of locations in a program by Loc = {loc 0 , loc 1 . . . , loc n-1 } where n is the number of allocation sites. In our representation, i.e., Figure 6, new locations are created by malloc operations. EXAMPLE 1. Figure 7 shows the control flow graph of the program seen in Figure 1. The two allocations at lines 17 and 18 are associated respectively with loc 0 and loc 1 .

Symbolic Range Analysis.

We start our pointer analysis by running an off-the-shelf range analysis parameterized on symbols. For the sake of completeness, we shall revisit the main notions associated with range analysis, which we borrow from Nazaré et al. [START_REF] Nazaré | Validation of memory accesses through symbolic analyses[END_REF]. We say that E is a symbolic expression, if and only if, E is defined by the grammar below. In this definition, s is a symbol and n ∈ N. The set of symbols s in a program form its symbolic kernel. The symbolic kernel is formed by names that cannot be represented as function of other names in the program text. Concretely, this set contains the names of global variables and variables assigned with values returned from library functions.

E ::= n | s | min(E, E) | max(E, E) | E -E | E + E | E/E | E mod E | E × E
We shall be performing arithmetic operations over the partially ordered set S = S E ∪ {-∞, +∞}, where S E is the set of symbolic expressions. The partial ordering is given by -

∞ < . . . < -2 < -1 < 0 < 1 < 2 < . . . + ∞.
There exists no ordering between two distinct elements of the symbolic kernel of a program. For instance, N < N +1 but there is no relationship between an expression containing N and another expression containing

M . A symbolic interval is a pair R = [l, u],
where l and u are symbolic expressions. We denote by R ↓ the lower bound l and R ↑ the upper bound u. We define the partially ordered set of (symbolic) intervals S 2 = (S × S, ), where the ordering operator is defined as:

[l 0 , u 0 ] [l 1 , u 1 ], if l 1 ≤ l 0 ∧ u 1 ≥ u 0
From the previous definitions, we define the semi-lattice SymbRanges of symbolic intervals as (S 2 , , , ∅, [-∞, +∞]), where the join operator " " is defined as:

[a 1 , a 2 ] [b 1 , b 2 ] = [min(a 1 , b 1 ), max(a 2 , b 2 )]
Our lattice has a least element ∅, such that:

∅ [l, u] = [l, u] ∅ = [l, u]
and a greatest element [-∞, +∞], such that:

[-∞, +∞] [l, u] = [l, u] [-∞, +∞] = [-∞, +∞]
For sake of clarity, we also define the intersection operator " ":

[a 1 , a 2 ] [b 1 , b 2 ] = ∅, if a 2 < b 1 or b 2 < a 1 [max(a 1 , b 1 ), min(a 2 , b 2 )], otherwise [-∞, +∞]
is absorbant and ∅ is neutral for . The result of range analysis is a function R : V → S 2 , that maps each integer variable i in a program to an interval [l, u], l ≤ u, e.g., R(i) = [l, u]. The more precise the technique we use to produce this result, the more precise our results will be. Nevertheless, the exact implementation of the range analysis is immaterial for the formalization that follows. In this paper, we are using the following widening operator on SymbRanges :

[l, u]∇[l , u ] =          [l, u] if l = l and u = u [l, +∞] if l = l and u > u [-∞, u] if l < l and u = u [-∞, +∞] if l < l and u > u
The only requirement that we impose on the implementation of range analysis is that it exists over SymbRanges , our lattice of symbolic intervals. We denote by (α SymbRanges , γ SymbRanges ) the underlying galois connection. EXAMPLE 2. A range analysis such as Nazaré et al.'s [START_REF] Nazaré | Validation of memory accesses through symbolic analyses[END_REF], if applied onto the program seen in Figure 3, will gives us that R(i

line 3 ) = [0, 0], R(i n.5 ) = [0, N -1], R(i n.7 ) = [0, N + 1].

Global Range Analysis of Pointers

As we have mentioned in Section 2 we use two different strategies to disambiguate pointers: the global and the local test. Our global pointer analysis goes over the entire code of the program, associating variables that have the pointer type with elements of an abstract domain that we will define soon. The local analysis, on the other hand, works only for small regions of the program text. We shall discuss the local test in Section 3.6. In this section, we focus on the global test, which is an abstract-interpretation based algorithm.

An Abstract Domain of Pointer Locations. We associate pointers with tuples of size n: (SymbRanges • ∪⊥) n ; n being the number of program sites where memory is allocated (the cardinal of Loc) and • ∪ is the disjoint union. Let @(loc i ) denotes the actual address value returned by the i th malloc of the program. By construction, all actual addresses are supposed to be offsets of a given @(loc i ). The abstract value GR(p) = (p 0 , . . . , p n-1 ) represents (an abstract version) of the set of memory locations that pointer variable p can address throughout the execution of a program: • p i = ⊥ if there is no adress in S which is an offset of @(loc i )

• p i = α SymbRanges ({d ∈ Z | s = @(loc i ) + d ∈ S}),
otherwise. The offsets from a given pointer are abstracted alltogether in the SymbRanges lattice.

The goal of our GR analysis is to compute such an abstract value for each pointer of the program. Some elements in a tuple GR(p) are bound to the undefined location, e.g., ⊥. These elements are not interesting to us, as they do not encode any useful information. Thus, to avoid keeping track of them, we rely on the concept of support, which we state in Definition 2. DEFINITION 2 (Support). We denote by supp GR (p) the set of indexes for which p i is not ⊥ :

supp GR (p) = {i | p i =⊥}.
For sake of readability, let us denote for instance GR(p) = (⊥, [l 1 , u 1 ], ⊥, [l 3 , u 3 ], ⊥), by the set GR(p) = {loc 1 + [l 1 , u 1 ], loc 3 + [l 3 , u 3 ]}. In the concrete world, this notation will mean that pointer p can address any memory location from @(loc 1 ) + l 1 to @(loc 1 ) + u 1 , and from @(loc 3 ) + l 3 to @(loc 3 ) + u 3 .

For instance, consider that l 1 = 3, u 1 = 5, l 3 = 3 and u 3 = 8. GR(p) = {loc 1 + [START_REF] Balakrishnan | Analyzing memory accesses in x86 executables[END_REF][START_REF] Bodik | ABCD: eliminating array bounds checks on demand[END_REF], loc 3 + [3, 8]} is then depicted in Figure 8. Now for the abstract operations: (⊥, . . . , ⊥) is the least element of our lattice, and

([-∞, ∞], . . . , [-∞, ∞]) the greatest one.
Given the two abstract values GR(p 1 ) = (p 1 0 , . . . p 1 n-1 ) and GR(p 2 ) = (p 2 0 , . . . p 2 n-1 ), the union GR(p 1 ) GR(p 2 ) is the tuple (q 0 , . . . , q n-1 ) where:

q i = ⊥ if p 1 i = p 2 i = ⊥ p 1 p 2 else
and GR(p1 ) GR(p 2 ) if an only if all involved (symbolic) intervals of p 1 are included in the ones of p 2 : ∀i ∈ [0..(n -1)], p 1 i p 2 i (considering ⊥ R and ⊥ R = R for all non-empty intervals R). We call this structure, formed by (SymbRanges • ∪ ⊥) n plus its partial ordering the lattice MemLocs . EXAMPLE 3. For the example depicted in Figure 7 where we only have two malloc sites denoted by loc 0 and loc 1 , we will obtain the following results: GR

(p) = GR(b) = {loc 0 + [0, 0]}, GR(m 0 ) = GR(s) = {loc 1 + [0, 0]}, GR(e) = loc 0 + [N, N ], GR(m 1 ) = loc 1 + [1, +∞], GR(i 7 ) = loc 0 + [N + strlen(m 0 ), N + strlen(m 0 ) + 1]
. How this mapping is found is discussed in the rest of this section.

Abstract semantics for GR, and concretisation. The abstract semantics of each instruction in our core language is given by Figure 9. Figure 9 defines a system of equations whose fixed point gives us an approximation on the locations that each pointer may dereference. We remind the reader of our notation: [l, u] ↓ = l, and [l, u] ↑ = u. In Figure 9, this notation surfaces in the semantics of intersections. The abstract interpretation of the pointer-related instructions in Figure 7 yields the results discussed in Example 3.

j : p = malloc (v) with v scalar ⇒ GR(p) = (⊥, . . . , [0, 0] j th component , . . .) p = free (v) with v scalar ⇒ GR(p) = (⊥, . . . , ⊥ ) v = v 1 ⇒ GR(v) = GR(v 1 ) q = p + c with c scalar variable ⇒     
GR(q) = (q 0 , . . . , q n-1 ) with

q i = ⊥ if p i = ⊥ p i + R(c) else q = φ(p 1 , p 2 ) ⇒ GR(q) = GR(p 1 ) GR(p 2 ) q = p 1 ∩ [-∞, p 2 ] ⇒     
GR(q) = (q 0 , . . . , q n-1 ) with

q i = ⊥ if (p 1 i = ⊥ or p 2 i =⊥) p 1 i [-∞, p 2 i ↑ ] else q = p 1 ∩ [p 2 , +∞] ⇒     
GR(q) = (q 0 , . . . , q n-1 ) with There remains to define how the abstract states will be concretised (@(loc i ) is the actual address returned by the i th malloc): DEFINITION 3 (Concretisation). Given GR(p) an abstract value (a set of "abstract addresses for p"), denoted by GR(p) = (p 0 , . . . , p n-1 ) then we define its concretisation as follows:

q i = ⊥ if (p 1 i = ⊥ or p 2 i = ⊥ ) p 1 i [p 2 i ↓ , +∞] else q = * p ⇒ GR(q) = ([-∞, ∞], . . . , [-∞, ∞]) * q = p ⇒ Nothing
γ(GR(p)) = i∈supp GR (p) {@(loc i ) + o, o ∈ p i }
The concretisation function of this abstract value is thus a set of (concrete) addresses, obtained by shifting a set of base addresses by a certain value in SymbRanges 1 . PROPOSITION 1. (α, γ) is a Galois connection. PROOF 1. Immediate since (α SymbRanges , γ SymbRanges ) is a galois connection. Solving the abstract system of contraints Following the abstract interpretation framework, we solve our system of constraints by computing for each pointer a growing set of abstract values until convergence.

However, as the underlying lattice SymbRanges has infinite height, widening is necessary to ensure that these sequence of iterations actually terminate. Our widening operation on pointers generalizes the widening operation on ranges. It is defined as follows: As usual, we only apply the widening operator on a cut set of the control flow graphs (here, only on φ functions).

Widening may lead our interpreter to produce very imprecise results. To recover part of this imprecision, we use a descending sequence of finite size: after convergence, we redo a step of symbolic evaluation of the program, starting from the value obtained after convergence. One example of analysis will be detailed later, in Section 3.9.

The abstract interpretation of loads and stores. In Figure 9, we chose not to track precisely the intervals associated with pointers stored in memory. In other words, when interpreting stores, e.g., q = * p, we assign the top value of our lattice to q. This decision is pragmatic. As we shall explain in Section 4, a typical compilation infra-structure already contains analyses that are able to track the propagation of pointer information throughout memory. Our goal is not to solve this problem. We want to deliver a fast analysis that is precise enough to handle C-style pointer arithmetic.

Answering GR Queries

Our queries are based on the following result, that is an immediate consequence of the fact our analysis is an abstract interpretation: PROPOSITION 2 (Correctness). Let p and p be two pointers in a given program then :

if supp GR (p) ∩ supp GR (p ) = ∅ or ∀i ∈ supp GR (p) ∩ supp GR (p ), p i p i = ∅ a 1 = malloc(2) a 2 = a 1 + 1 a 3 = ɸ(a 1 , a 2 ) a 4 = a 3 + 1 a 5 = a 3 + 2 GR(a 1 ) = {loc 1 + [0, 0]} GR(a 2 ) = {loc 1 + [1, 1]} GR(a 3 ) = {loc 1 + [0, 1]} GR(a 4 ) = {loc 1 + [1, 2]} GR(a 5 ) = {loc 1 + [2, 3]} LR(a 1 ) = loc 1 + [0, 0] LR(a 2 ) = loc 1 + [1, 1] LR(a 3 ) = loc 2 + [0, 0] LR(a 4 ) = loc 2 + [1, 1] LR(a 5 ) = loc 2 + [2, 2]

Global Analysis

Local Analsys Figure 10. Example that illustrates imprecision of global analysis due to lack of path-sensitiveness.

then γ(GR(p)) ∩ γ(GR(p )) = ∅.
In other words, if the abstract values of two different pointers of the program have a null intersection, then the two concrete pointers do not alias. This result is directly implied by the abstract interpretation framework. Thanks to this result, we implement the query Q GR (p, p ) as :

• If GR(p) and GR(p ) have an empty intersection, then "they do not alias". • Else "they may alias".

Local Range Analysis of Pointers

The global pointer analysis is not path sensitive. As a consequence, this analysis cannot, for instance, distinguish the effects of different iterations of a loop upon the actual value of a pointer, or the effects of different branches of a conditional test on that very pointer. The program in Figure 10 illustrates this issue. Pointers a 4 and a 5 clearly must not alias. Yet, their abstract states have non-empty intersections for loc 1 . Therefore, the query mechanism of Section 3.5 would return a "may-alias" result in this case.

To solve this problem, we have developed a local version of our pointer analysis. We call it local because it creates new locations for every φ-function. Our local range analysis is simpler than its global counterpart. We solve it in a single iteration of abstract interpretation applied on the instructions of our core language. Instructions are evaluated abstractly in the order given by the program's dominance tree. Figure 11 gives the abstract semantics of each instruction. The abstract value LR(p) exists in (Loc∪N ewLocs)× SymbRanges where N ewlocs denotes a set of "fresh location variables", that are computed by invocation of the function N ewLocs(). As before, we write loc + R instead of (loc, R). Similarly to γ GR , γ LR denotes the set of abstract addresses from @(loc) + R ↓ to @(loc) + R ↑ .

To find a solution to the local analysis, we solve the system provided by the abstract rules seen in Figure 11. This resolution process involves computing an increasing sequence of abstract values for each pointer p of the program. Contrary to the global analysis, this analysis is based on a finite lattice, we do not need any widening operator. Figure 10 (Right) shows the result of the local analysis. Contrary to the

p = malloc (v) with v scalar ⇒ LR(p) = N ewLocs() + [0, 0] v = v 1 ⇒ LR(v) = LR(v 1 ) q = p + c
with c scalar variable and

LR(q) = loc i + [l, u] ⇒ LR(q) = loc i + ([l, u] + R(c)) j : q = φ(p 1 , p 2 ) ⇒ LR(q) = N ewLocs() + [0, 0] q = p 1 ∩ [-∞, p 2 ] q = p 1 ∩ [p 2 , +∞] ⇒ LR(q) = LR(p 1 ) q = * p 1 ⇒ LR(q) = N ewLocs() + [0, 0] * q = p 1 ⇒

Answering LR Queries

The correction for the local analysis is stated by the following proposition: Thanks to this result, we implement the query Q LR (p, p ):

• If LR(p) and LR(p ) have a common base pointer with ranges that do not intersect, then "they do not alias". • Else "they may alias".

Complexity

The e-SSA representation ensures that we can implement our analysis sparsely. Sparsity is possible because the e-SSA form renames variables at each program point where new abstract information, e.g., ranges of integers and pointers, arises. According to Tavares et al. [START_REF] Tavares | Parameterized construction of program representations for sparse dataflow analyses[END_REF], this property -single information -is sufficient to enable sparse implementation of non-relational static analyses [START_REF] Tavares | Parameterized construction of program representations for sparse dataflow analyses[END_REF]. Therefore, the abstract state of each variable is invariant along the entire live range of that variable. Consequently, the space complexity of our static analysis is O(|V | × I), where V is the set of names of variables in the program in e-SSA form, and I is a measure of the size of the information that can be bound to each variable.

We apply widening after one iteration of abstract interpretation. Thus, we let the state of a variable to change first |). This observation also prevents our algorithm from generating expressions with very long chains of "min" and "max" expressions. Therefore, I, the amount of information associate with a variable, can be represented with O(1) space. As a consequence of this frugality, our static analysis runs in O(|V |) time, and requires O(|V |) space.

A wrap-up Example

Example 4 shows how our analysis works on the program seen in Figure 1. EXAMPLE 4. Figure 7 shows the control flow graph (CFG) of the program in Figure 1. Our graph is in e-SSA form [START_REF] Bodik | ABCD: eliminating array bounds checks on demand[END_REF]. Figure 12 shows the result of widening ranges after one round of abstract interpretation (stabilization achieved), and a descending sequence of size two. Our system stabilizes after each instruction is visited four times. The first visit does initialization, the second widening (and stabilization check), and the last two build the descending sequence.

This example illustrates the need of widening to ensure termination. Our program has a cycle of dependencies between pointers i 1 , i 2 and i 3 . If not for widening, pointer i 3 , incremented in line 5 of Figure 1 would grow forever. Thus, as in Abstract Interpretation, we must break the cyclic dependences between our pointers under analysis, by means of insertion of widening points (identify points in the CFG where to apply widening to insure convergence).

Returning to our example of Figure 1, we are interested in knowing, for instance, that the memory access at line 6 is independent on the accesses that happen at line 10. To achieve this goal, we must bound the memory regions covered by pointers i 3 and i 7 . A cyclic dependence happens at the operation i++, because in this case, we have a pointer being used as both, source and destination of the update. Thus, we should have inserted a widening point at stores and load instructions. However, in the Abstract Interpreter depicted in Figure 9, it was sufficient to insert widening points at φ functions (as we already said before) because :

• heads of loops are φ functions (thus dependencies between variables of different iteration of loops are broken). • we are working on (e-)SSA form programs; thus, the only inter-loop dependencies are successive stores to the same variable : *q=..., *q=.... The value GR(q) is the union of all information gathered inside the loop. (In essence, memory addresses are not in static single assignment form, i.e., we could have the same address being used as the target of a store multiple times). This information might grow forever; hence, we would have inserted a widening point on the last write. In our case, the information we store is already the top of our lattice; hence, there is no need for widening.

Var GR LR b, p, i 0 ([0, 0], ⊥) loc 0 + [0, 0] m 0 , s (⊥, [0, 0]) loc 1 + [0, 0] i 1 ([0, 0], ⊥) loc 2 + [0, 0] i 2 ([0, 0], ⊥) loc 2 + [0, 0] t 0 ([1, 1], ⊥) loc 2 + [1, 1] e ([N, N ], ⊥) loc 0 + [N, N ] Starting i 3 ([1, 1], ⊥) loc 2 + [1, 1] state i 4 (⊥, ⊥) loc 2 + [0, 0] f ([k, k], ⊥) loc 0 + [k, k] m 1 (⊥, [0, 0]) loc 3 + [0, 0] m 2 (⊥, [1, 1]) loc 3 + [1, 1] i 5 (⊥, ⊥) loc 4 + [0, 0] i 6 (⊥, ⊥) loc 4 + [0, 0] i 7 (⊥, ⊥) loc 4 + [1, 1] i 1 ([0, +∞], ⊥) i 2 ([0, +∞], ⊥) t 0 ([1, +∞], ⊥) Growing i 3 ([1, +∞], ⊥) iterations i 4 ([N, +∞], ⊥) + widening m 1 (⊥, [0, +∞]) m 2 (⊥, [1, +∞]) i 5 ([N, +∞], ⊥) i 6 ([N, k -1], ⊥) i 7 ([N + 1, k], ⊥) i 2 ([0, N -1], ⊥) t 0 ([1, N ], ⊥) After one i 3 ([1, N ], ⊥) descending m 1 (⊥, [0, +∞]) step m 2 (⊥, [1, +∞]) i 1 ([0, N ], ⊥) After two i 3 ([1, N ], ⊥) descending i 4 ([N, N ], ⊥) steps i 5 ([N, k], ⊥) i 6 ([k -1, k], ⊥) i 7 ([k, k + 1], ⊥) Figure 12
. Abstract interpretation of CFG seen in Figure 7 (program in Figure 1). For GR, we associate loc 0 with the malloc at line 17 and loc 1 with the malloc at line 18 (of the program). Only changes in GR and LR are rewritten after the growing and descending iterations. We let k = N+strlen(m 0 ).

Experiments

We have implemented our range analysis in the LLVM compiler, version 3.5. In this section, we show a few numbers that we have obtained with this implementation. All our experiments have been performed on an Intel i7-4770K, with 8GB of memory, running Ubuntu 14.04.2. Our goal with these experiments is to show: (i) that our alias analysis is more precise than other alternatives of practical runtime; and (ii) that it scales up to large programs.

On the Precision of our Analysis. In this section, we compare our analysis against the other pointer analyses that are available in LLVM 3.5, namely basic and SCEV. The first of them, although called "basic", is currently the most effective alias analysis in LLVM, and is the default choice at the -O3 optimization level. It relies on a number of heuristics to disambiguate pointers2 :

• Distinct globals, stack allocations, and heap allocations can never alias. • Globals, stack allocations, and heap allocations never alias the null pointer. • Different fields of a structure do not alias.

• Indexes into arrays with statically differing subscripts cannot alias. • Many common standard C library functions never access memory or only read memory. • Function calls cannot reference stack allocations which never escape from the function that allocates them As we see from the above list, the basic alias analysis has some of the capabilities of the technique that we present in this paper, namely the ability to distinguish fields and indices within aggregate types. In this case, such disambiguation is only possible when the aggregates are indexed with constants known at compilation time. For situations when these indices are symbols, LLVM relies on a second kind of analysis to perform the disambiguation: the "scalar-evolutionbased" (SCEV) alias analysis. This analysis tries to infer closed-form expressions to the induction variables used in loops. For each loop such as:

for (i = B; i < N; i += S) { ... a[i] ... }
this analysis associates variable i with the expression i = B+iter ×S, i ≤ N. The parameter iter represents the current iteration of the loop. With this information, SCEV can track the ranges of indices which dereference array a within the loop. Contrary to our analysis, SCEV is only effective to disambiguate pointers accessed within loops and indexed by variables in the expected closed-form.

Figure 13 shows how the three different analyses fare when applied on larger benchmarks. For this experiment we have chosen three benchmarks that have been used in previous work that compares pointer analyses: Prolangs [START_REF] Ryder | A schema for interprocedural modification side-effect analysis with pointer aliasing[END_REF], PtrDist [START_REF] Zhao | Dynamic memory optimization using pool allocation and prefetching[END_REF] and MallocBench [START_REF] Grunwald | Improving the cache locality of memory allocation[END_REF]. We first notice that in general all the pointer analyses in LLVM disambiguate a relatively low number of pointers. This happens because many pointers are passed as arguments of functions, and, not knowing if these functions will be called from outside the program, the analyses must, conservatively, assume that these parameters may alias. Second, we notice that our pointer analysis is one order of magnitude more precise than the scalar-evolution based implementation available in LLVM. Finally, we notice that we are able to disambiguate more queries than the basic analysis. Furthermore, our results complements it in non-trivial ways. In total, we tried to disambiguate 3.093 million pairs of pointers. Our analysis found out that 1.29 million pairs reference non-overlapping regions. The basic analysis has been able to distinguish 953 thousand pairs. By combining these two analyses, we extended this number to 1.439 thousand pairs of pointers. SCEV could not increase this number any further.

Figure 14 shows the proportion of queries that we have been able to disambiguate with the global test of Section 3.4. On the Scalability of our Analysis. The chart in Figure 15 shows how our analysis scales when applied on programs of different sizes. We have used the 50 largest programs in the LLVM benchmark suite. These programs gave us a total of 800,720 instructions in the LLVM intermediate representation, and a total of 241,658 different pointer variables. We analyzed all these 50 programs in 8.36 seconds. We caneffectively -analyze 100,000 instructions in about one second. In this case, we are counting only the time to map variables to values in SymbRanges . We do not count the time to query each pair of pointers, because usually compiler optimizations perform these queries selectively, for instance, only for pairs of pointers within a loop. Also, we do not count the time to run the out-of-the-box implementation of range analysis mentioned in Section 3.3, because our version of it is not implemented within LLVM. It runs only once, and we query it afterwards, never having to re-execute it. The chart in Figure 15 provides strong visual indication of the linear behavior of our algorithm. We have found, indeed, cogent evidence pointing in this direction. The linear correlation coefficient (R) indicates how strong is a linear relationship between two variables. The closer to one, the more linear is the correlation. The linear correlation between time and number of instructions for the programs seen in Figure 15 is 0.982, and the correlation between time and number of pointers is 0.975.

Related Work

The contribution of this work is a new representation of pointers, based on the SymbRanges lattice, and an algorithm to reach a fixed point in this lattice, based on abstract in-terpretation. This contribution complements classic work on pointer analysis. In other words, our representation of pointers can be used to enhance the precision of algorithms such as Steensgard's [START_REF] Steensgaard | Points-to analysis in almost linear time[END_REF], Andersen's [START_REF] Andersen | Program Analysis and Specialization for the C Programming Language[END_REF], or even the state-ofthe-art technique of Hardekopf and Lin [START_REF] Hardekopf | Flow-sensitive pointer analysis for millions of lines of code[END_REF]. These techniques map pointers to sets of locations, but they could be augmented to map pointers to sets of locations plus ranges. Furthermore, the use of our approach does not prevent the employment of acceleration techniques such as lazy cycle detection [START_REF] Hardekopf | The ant and the grasshopper: fast and accurate pointer analysis for millions of lines of code[END_REF], or wave propagation [START_REF] Pereira | Wave propagation and deep propagation for pointer analysis[END_REF].

There exist previous work that used similar lattices as ours, albeit different resolution algorithms. For instance, much of the work on automatic parallelization has some way to associate symbolic offsets, usually loop bounds, with pointers. Michael Wolfe [START_REF] Wolfe | High Performance Compilers for Parallel Computing[END_REF]Ch.7] and Aho et al. [START_REF] Aho | Compilers: Principles, Techniques, and Tools[END_REF]Ch.11] have entire chapters devoted to this issue. The key difference between our work and this line of research is the algorithm to solve pointer relations: they resort to integer linear programming (ILP) or the Greatest Common Divisor test to solve diophantine equations, whereas we do abstract interpretation. Even Rugina and Rinard [START_REF] Rugina | Symbolic bounds analysis of pointers, array indices, and accessed memory regions[END_REF], who we believe is the state-of-the-art approach in the field today, use integer linear programming to solve symbolic relations between variables. We speculate that the ILP approach is too expensive to be used in large programs; hence, concessions must be made for the sake of speed. For instance, whereas the previous literature that we know restrict their work to pointers within loops, we can analyze programs with over one million assembly instructions in a few seconds.

There exist work that, like ours, also associates intervals with pointers, and solves static analysis via abstract interpretation techniques. However, to the best of our knowledge, these approaches have a fundamental difference to our work: they use integer intervals a là Cousot [START_REF] Cousot | Abstract interpretation: a unified lattice model for static analysis of programs by construction or approximation of fixpoints[END_REF], whereas we use symbolic intervals. The inspiration for much of this work springs from Balakrishnan and Reps notion of Value Set Analysis [START_REF] Balakrishnan | Analyzing memory accesses in x86 executables[END_REF]. Integer intervals have also being used by Yong et al. [START_REF] Yong | Pointer-range analysis[END_REF] and, more recently, by Oh et al. [START_REF] Oh | Selective contextsensitivity guided by impact pre-analysis[END_REF]. In the latter case, Oh et al. use pointer disambiguation incidentally, to demonstrate their ability to implement efficiently static analyses in a context-sensitive way. Even though integer ranges fit well the need of machine code, as demonstrated by Balakrishnan and Reps, we believe that further precision requires more expressive lattices. We have not implemented value set analysis, but we have tried a simple experiment: we counted the number of pointers that have integer ranges, and compared this number against the quantity of pointers that have symbolic ranges. We found out that 20.47% of the pointers in our three benchmark suites have exclusively symbolic ranges. Classic range analysis would not be able to distinguish them. Notice that numeric ranges are more common among pointer variables than among integer, because fields within structs -a very common construct in Care indexed through integers. Finally, the fact that we use

Figure 1 .

 1 Figure 1. Example of program that builds up messages as sequences of serialized bytes. We are interested in disambiguating the locations accessed at lines 6 and 10.

Figure 2 . 8 / 9 Figure 3 . 10 Figure 4 .

 2893104 Figure 2. Array p in the routine prepare seen in Fig.1. Lines 6 and 10 represent the different stores in the figure.

Figure 5 .

 5 Figure 5. Overview of our pointer analysis.

Figure 7 .

 7 Figure 7. Control flow graph of program seen in Figure 1

Figure 8 .

 8 Figure 8. The concrete semantics of GR(p) = {loc 1 + [3, 5], loc 3 + [3, 8]}. Dark grey cells denote possible (concrete) values of p.

Figure 9 .

 9 Figure 9. Constraint generation for GR with GR(p) = (p 0 , . . . , p n-1 ) given p in the right hand side of rules

DEFINITION 4 .

 4 Given GR(p) and GR(p ) with GR(p) GR(p ), we define the widening operator: GR(p)∇GR(p ) = (p 0 ∇p 0 , . . . , p n-1 ∇p n-1 ), where ∇ denotes the widening on SymbRanges , extended with ⊥∇⊥ = ⊥ and ⊥∇[l, u] = [l, u].

PROPOSITION 3 (

 3 Correctness). Let p and p be two pointers in a given program, and γ LR be the concretization of the abstract map LR, which we state like in Definition 3. If LR(p) = loc + R and LR(p ) = loc + R , then if loc = loc and R R = ∅ then γ(LR(p))∩γ(LR(p ) = ∅. In other words, p and p never alias.

Author

  Version of « Symbolic Range Analysis of Pointers » Published in Code Generation and Optimization, Barcelona, 2016 from [⊥, ⊥] to [s l , s u ], where s l = -∞, and s u = +∞. From there, we can reach either [-∞, s u ] or [s l , +∞]. And, finally, this abstract state can jump to [-∞, +∞]. Hence, our time complexity is O(3 × |V |) = O(|V

  Nothing Figure 11. Constraint generation for LR global analysis, we have a new location bound to variable a 3 , which is defined by a φ operator. The range of this new location is [0, 0]. The other variables that are functions of a 3 , e.g., a 4 and a 5 , have now non-intersection ranges associated with this new memory name.

  Figure[START_REF] Nazaré | Validation of memory accesses through symbolic analyses[END_REF]. Runtime of our analysis for the 50 largest benchmarks in the LLVM test suite. Each point on the X-axis represents a different benchmark. Benchmarks are ordered by size. This experiment took less than 10 seconds.ing the execution of the program. The local test has been able to disambiguate 6.55% of all the addresses used in our benchmarks. The rest of the disambiguation was obtained by comparing off-sets from different locations.

	1.E+05
	1.E+04
	1.E+03
	1.E+02
	1.E+01
	0	5	10	15	20	25	30	35	40	45	50
		#Instruc3ons	#Pointers		Run3me (ms)

While speaking about symbolic ranges, we also have to concretize the values involved in the bounds of p i , that is we shall use the actual values between S(p i↓ ) and S(p i↑ ).

This list has been taken from the LLVM documentation, available at http://llvm.org/docs/AliasAnalysis.html in September of 2015

Acknowledgment: This project is supported by the Brazilian Ministry of Science and Technology through CNPq, the Intel Corporation (the ISRA eCoSoC project) and the INRIA-FAPEMIG cooperation grant (The Prospiel project).

* This work was partially supported by the LABEX MILYON (ANR-10-LABX-0070) of Université de Lyon, within the program "Investissement d'Avenir" (ANR-11-IDEX-0007) operated by the French National Research Agency (ANR) Categories and Subject Descriptors D -Software [D.3

Author Version of « Symbolic Range Analysis of Pointers » Published in Code Generation and Optimization, Barcelona, 2016

Bodik's e-SSA form [START_REF] Bodik | ABCD: eliminating array bounds checks on demand[END_REF] distinguish our abstract interpretation algorithm from previous work. This representation lets us solve our analysis sparsely, whereas Balakrishnan's algorithm works on a dense representation that associates facts with pairs formed by variables and program points.

Conclusion

In this paper we have presented a new alias analysis technique that handles, within the same theoretical framework, the subtleties of pointer arithmetic and memory indexation. Our technique can disambiguate regions within arrays and C-like structs using the same abstract interpreter. We have achieved precision in our algorithm by combining alias analysis with classic range analysis on the symbolic domain. Our analysis is fast, and handles cases that the implementations of pointer analyses currently available in LLVM cannot deal with. In future work, we plan to investiguate better splitting strategies and other more expressive lattices to improve the global precision of our analyses.


Yasser Moursy
email: yasser.moursy@lip6.fr

Hao Zou

Ramy Iskander

Pierre Tisserand

Dieu-My Ton

Giuseppe Pasetti

Ehrenfried Seebacher

Alexander Steinmair

Thomas Gneiting

Heidrun Alius

Towards Automatic Diagnosis of Minority Carriers Propagation Problems in HV/HT Automotive Smart Power ICs

Keywords: smart power integrated circuits, substrate modeling, substrate noise, minority carriers I

In this paper, a proposed methodology to identify the substrate coupling effects in smart power integrated circuits is presented. This methodology is based on a tool called AUTOMICS to extract substrate parasitic network. This network comprises diodes and resistors that are able to maintain the continuity of minority carrier concentration. The contribution of minority carriers in the substrate noise is significant in high-voltage and high temperature applications. The proposed methodology along with conventional latch-up problem identification for a test case automotive chip AUTOCHIP1 are presented. The time of the proposed methodology is significantly shorter than the conventional one. The proposed methodology could significantly shorten the time-to-market and ameliorate the robustness of the design.

INTRODUCTION

Automotive applications have become a growing market for smart power integrated circuits. Nowadays, the automotive industry aims at developing green cars along with ameliorating network connectivity. Hence, the smart power ICs are becoming increasingly attractive to miniaturize the systems and implement more functions in the vehicles. Consequently, addressing electronics failures has a spate of interest in order to guarantee more robust and reliable products. Currently, the design centers can prove by simulation the functional modes. As a result, the number of the electronic failures are drastically reduced thanks to Design for Reliability (DfR) method [START_REF] Ognjanovic | Design for reliability based methodology for automotive gearbox load capacity identification[END_REF] [START_REF] Sasun | Physics of failure based reliability prediction method for two-wheeler[END_REF] globally deployed in automotive besides, the monolithic integration in ICs where reliability is ameliorated by manufacturing.

However, failures due to substrate noise coupling are still reported in tests after fabrication. These failures could be latchup that can lead to the destruction of the circuit. The designer can extract them manually and consequently, the robustness of the circuit depends on the background of the designer. In [START_REF] Schenkel | Substrate current effects in smart power ICs[END_REF], the noise in smart power ICs is addressed and a TCAD-guided design methodology for substrate current control is proposed.

Addressing more complex structure with this methodology could be impractical as it takes excessively long time. In [START_REF] Conte | A circuit-level substrate current model for smart-power ICs[END_REF], a modeling methodology for the PN junctions contributing in the substrate noise has been introduced. Additional terminals to the conventional diodes and resistors are used to preserve the minority carrier diffusion. The key advantage of this technique is its compatibility with circuit simulator. Hence, the simulation time is considerably shorter than TCAD simulations. However, the substrate network is extracted manually that is tedious for complex circuit analyses.

The main objective of this paper is to propose an automated approach for modeling the substrate. This approach relies on automatically extracting the substrate parasitic components from the layout. These components are modeled using special types of diodes and resistors that ensure the continuity of the minority carrier propagation [START_REF] Stefanucci | Modeling Minority Carriers Related Capacitive Effects for Transient Substrate Currents in Smart Power ICs[END_REF] [6] [START_REF] Stefanucci | Substrate modeling to improve reliability of high voltage technologies[END_REF]. The tool used for the extraction and composing the equivalent substrate network is called AUTOMICS. We illustrate the conventional method for failure detection in the industrial verification flow and how our proposed verification flow could shorten the time-to-market and ameliorate the design.

The paper is organized as follows. In Section II, the AUTOMICS tool flow is illustrated and the calibration of the tool for the parameters of the process technology. The proposed design/verification flow is explained in Section III. Section IV shows a study for an industrial test case chip, finally, the conclusions are drawn in Section V.

II. SUBSTRATE PARASITIC EXTRACTION WITH AUTOMICS TOOL

AUTOMICS tool is capable of modeling the substrate considering the minority and majority carriers propagation. It depends on special types of diodes and resistors presented in [START_REF] Stefanucci | Modeling Minority Carriers Related Capacitive Effects for Transient Substrate Currents in Smart Power ICs[END_REF]. The extraction procedure has two steps:

A. Geometrical features extraction using AUTOMICS tool, B. Technological parameters calibration.

A. AUTOMICS tool geometrical features extraction

The input to AUTOMICS tool is the layout of the design as shown in Figure 1. As layout designs are in two dimensions (xand y-axes), the z-axis is technology dependent. For instance, the depth of mask layers, and their doping profiles. Thus, for extraction, this information is given as guide rules for extraction. These extraction rules are specified for each process technology. The entire Computer-Aided-Design (CAD) framework is developed using the "OpenAccess" open source API [START_REF] Zou | A novel CAD framework for substrate modeling[END_REF].

The processes operated on the layout are reduction, meshing, and extraction. In the reduction phase, a reduction mechanism is used to filter out the layers that are not contributing in the substrate parasitic coupling, e.g. metal layers. The meshing phase, a 2D meshing strategy is applied to the reduced layout. In conjunction with the extraction rules, a substrate 3D mesh is constructed. A non-uniform meshing mechanism is used to reduce the complexity. The extraction phase, parasitic devices are extracted and classified based on the doping profiles and orientation [START_REF] Zou | Substrate noise modeling with dedicated CAD framework for smart power ICs[END_REF]. A netlist describing the extracted parasitic components is generated. In the postprocessing stage, the extracted substrate netlist is ready for simulation with or without original schematic netlist. Simulations with the original netlist could be used to visualize the effect of minority carriers propagation in the substrate.

B. Model parameters calibration

Model calibration is a necessary step to estimate the doping concentrations and generate a parasitic model library for a certain technology. In general, TCAD simulations can be used in such calibration process. However, this approach needs a detailed input with respect to doping concentrations and doping profiles, which may not be convenient in industrial test cases. Therefore, a special set of test structures was implemented to characterize and model the basic diode behavior as well as the interaction between different well implants due to minority carriers taking into account shape and distance parameters of the well dimensions. The simulations are compared to measured data and an optimization algorithm is applied to tune the available model parameters to the measured curves. The calibration process in principle is shown in Figure 2.

III. PROPOSED SMART POWER ICS DESIGN METHODOLOGY

Certain design failures due to minority carriers coupling are only discovered after fabrication. The conventional design flow is shown in Figure 3, after fabrication an evaluation and validation step is done. Since the automotive applications are safety critical, extensive number of tests are performed on the chip. During these harsh tests, the design may fail. For welldesigned circuits, the simulator shows no errors which may not be the real case. Thus, parasitic hand extractions should be done considering the parasitic components that normally are ignored by the current commercial tools e.g. lateral NPN transistors, that is a layout-dependent, can cause a latch-up failure. After identifying the source of failure, the circuit is redesigned and refabricated. In some cases, this loop could be repeated more than one time, prolonging the time-to-market.

Our proposal is to add extra steps to the conventional flow before fabrication. In these steps, the substrate parasitic components are extracted from the layout using AUTOMICS tool. The output of the tool is substrate equivalent network that considers the majority and minority carriers propagation. Then, it can be simulated with the original design to visualize the effects of minority carriers propagation in the substrate. DC analyses can be done to estimate the amount of dc leakage current coupled to different wells on the substrate. This estimation can be used to reduce coupling between different wells by inserting guard rings or increase the distance. In addition to DC analyses, transient simulations could show effects due to junction capacitances of different wells on the same substrate. These capacitances significantly affect the system operation when abrupt voltage changes occur. Such effects could be the main source of latch-up in some circuits as will be discussed in next section.

The proposed methodology has the potential to reduce the number of redesigns. As the designer would be able to check the design robustness against effects that do not appear in the circuit simulations only. Furthermore, he would be able to optimize the isolation techniques e.g. guard rings, hence reducing the silicon area. These two factors can lead to reduce the cost and shorten the time-to-market. The key advantage is the short simulation time along with acceptable accuracy. The AUTOCHIP1 is a smart power-IC developed in 0.35um HV-CMOS technology from AMS. It has been developed for 24V Valeo automotive technology applications and is able to sustain overvoltages up to 90V. It aims at controlling 28V up to 32V heavy duty alternators by rotor high current driver up to 5A in order to provide alternator output current to the battery and loads. It integrates on the same die a DC-DC converter to drive rotor coil, similar to the one described in [START_REF] Saponara | HV-CMOS Design and Characterization of a Smart Rotor Coil Driver for Automotive Alternators[END_REF], different I/O interfaces, a power management block and a small digital part that manage a serial communication with an external digital core, implemented on an FPGA. The FPGA solution has been chosen to increase the flexibility, like in [START_REF] Fanucci | An high voltage CMOS voltage regulator for automotive alternators with programmable functionalities and full reverse polarity capability[END_REF].

This paper focuses on studying the DC-DC converter block. It includes a high-side switch and a low-side freewheeling diode. This stage regulates the current of the alternator rotor. One important point is the electrical connection and disconnection of the product that has high dV/dt up to 20V/ns that can be intermittent or during maintenance. In this section, the conventional procedure used to identify a problem in AUTOCHIP1 and our proposed methodology is discussed.

A. Conventional problem identificaion

The complete evaluation of the AUTOCHIP1 performances has been done in AMS and VALEO laboratories in parallel, performing the same tests and then cross-checking the results. One of the tests was the fast battery connection. Many fast battery connections at different temperatures and different voltages have been done. Through these tests, AUTOCHIP1 was very sensitive. As soon the battery is connected to a supply voltage of 50V, the device was severely damaged. Figure 4 shows the microphotograph of the AUTOCHIP1 after the test. The origin of this issue was a latch-up of the supply connection due to the activation of a parasitic SCR structure. Long investigations have taken more than three months, involving also emission microscope EMMI to localize and identify the pn junctions in forward direction during the test, and substrate parasitic estimation. A new schematic has been done with involved substrate parasitics added by hand and simulated in order to replicate the measurement behavior. Since hand-made parasitic extraction was needed, as shown in Figure 5, many iterations were needed to rebuild the entire involved substrate network and to replicate correctly the measurement behavior. Such simulations are needed to increase the confidence level of the redesign and to validate the root cause explanation.

After reasoning and validation, a re-design to avoid the triggering of the parasitic SCR structure. AUTOCHIP2 was designed and simulated in conjunction with the hand-extracted substrate network. The fabricated AUTOCHIP2 manage to pass the battery connection test safely.

B. Proposed problem identification

Our proposed methodology can be used to identify such problem. At this stage, the layout of the design is available and can be processed by AUTOMICS tool. First, the calibration for the technology parameters is done to estimate the doping concentrations. Then, a reduction rules file is written to consider only the layers contributing in the substrate noise. For simplicity, only the deep Nwells and substrate contacts of the PMOS transistor are considered. The tool output is a schematic of the substrate parasitics. This schematic has three main terminals bulk of the PMOS, the Ntype guard ring of the PMOS and the p-substrate contact. Figure 6 shows the testbench used for simulating AUTOCHIP1. It is noted that there is no connections for the guard rings in the design schematic as it is considered as a layout feature. The equivalent resistances of the contacts are estimated and inserted in the schematic. The voltage source V1 is set to 0V then ramps up to 12V in 10ns. As depicted in Figure 7, in case of circuit schematic only, a current spike appears at the supply connection and then the current returns to zero in fraction of microsecond. This indicates that the chip works properly. However, when the extracted substrate network is connected and the simulation is repeated. A holding current appears indicating the presence of a latched SCR structure. The results show an acceptable agreement with the measurement on the AUTOCHIP1. The simulation time and number of extracted components are shown in Table 1. The simulations were performed with a CPU sixcore AMD Opteron (tm) Processor.

These results could be used by a designer to analyze the latchup structure and provide solution for this problem before fabrication. Obviously, using such tool would save long time for fabrication and redesign and consequently lowering the cost and ameliorate the circuit design.

V. CONCLUSION

In this paper, we proposed a new methodology for analyzing smart power ICs. This methodology is based on extracting substrate parasitic components using a tool called AUTOMICS. The substrate is modeled as a network of diodes and resistors that preserve the continuity of the minority carrier propagation in the substrate. This aspect is essential for substrate analysis especially in case of high voltage and high temperature circuits. A test case of a chip called AUTOCHIP1 is presented. The parasitic hand-extraction took months plus the redesign and fabrication time. AUTOCHIP1 has been studied by the tool and the latch-up was detected. The proposed methodology gives an insight of the substrate behavior and shorten the time to market and hence reduces the cost.

Figure 1

 1 Figure 1 proposed substrate parasitic extraction flow

  Figure 3 proposed design flow with AUTOMICS

Figure 6

 6 Figure 6 AUTOCHIP1 simulation testbench

Table 1

 1

	substrate simulation parameters
	Number of extracted components	15083
	Simulation CPU time (s)	12.5
	Simulation elapsed time (s)	30.7

This work was supported by the European commission through the European FP7 within the AUTOMICS Project.


-1-

NEW APPROACH TOSCALING RULES FOR SOLAR AND PLANETARY DYNAMOS

Bertrand C. Barrois 1 Glorified dimensional analysis is used to derive scaling rules for internal and external magnetic field strengths and various time scales. Naive dimensional analysis is inconclusive because of multiple time scales, but physical arguments serve to weed out irrelevant parameters. Time scales can be derived from linearized instability analysis instead of ill-founded assumptions of Magnetic-Archimedean-Coriolis (MAC) balance. Further relationships can be derived from high-level models of coupled main field components and differential rotation. The ratios of the external dipole field to internal magnetic fields and of differential to overall rotation depend on details of the dynamo mechanism.

The dynamo mechanism by which turbulent convection sustains the magnetic fields of stars and planets has never been explained in satisfactory detail, but it is slowly yielding to numerical simulation. The goal of the present paper is to derivepredictive scaling rules for internal and external magnetic field strengths and important time scales in terms of physically relevant parameters.

Numerous authors (e.g., [START_REF] Christensen | Dynamo scaling laws and applications to the planets[END_REF][START_REF] Christensen | Energy flux determines magnetic field strength of planets and stars[END_REF]Christensen , 2010;;[START_REF] Davidson | Scaling laws for planetary dynamos[END_REF]Oruba 2014) have addressed the scaling problem by various theoretical and numerical methods, but no consensus has emerged. There appear to be significant differences between stars and planets, as well as between fast and slow rotation regimes.

Naive dimensional analysis is inconclusive because there are no less thanfourdifferentcombinations with units of time:

 Rotation period of the body  Inductive decay time for the magnetic dipole field (due to diffusion)  Mechanical time scale derived from the luminosity (or other power source)  Magnetic time scale derived from the magnetic energy density We shall seek a scaling law of the general form: ???

??? ??? One possible clue to differences between the Earth and Sun is the ratio / Q TT  -- roughly 400 for the Earth, but only 0.7 for the Sun. Other possibly significant differences involve the dissipation number and the thickness of the convective shell relative to radius.

~( / ) BQ T T T T f   
Besides time scales, there are dimensionless parameters of possible relevance:  The "dissipation number" for convective heat transport, which says how many times heat supplied at the base of the convective shell can be converted to kinetic energy by buoyancy and then redissipated on its way out.  The ratio of inner to outer radii of the convective shell  Ratios among diffusion coefficients (e.g., Prandtl number)

The lore of turbulent phenomena teaches that at extreme Reynolds numbers, kinematic viscosity (i.e., velocity diffusion) is virtually irrelevant in the inertial regime. By the same token, magnetic diffusion ( 1/ ) B    ought to be irrelevant if the time -3-scale for inductive decay is very long. The electrical conductivity seems relevant only to boundary conditions.

Within the Sun, the electrical conductivity of plasma scales as 3/2 () Temperature , according to the Spitzer-Härm relationship (Spitzer 1954[START_REF] Mitchner | Predictive scaling laws for rotating spherical dynamos[END_REF].The diffusion time is so long that time-dependent fields cannot penetrate the non-convective inner zone.

In the Earth, the conductivity of molten iron mixtures is uncertain, but it seems that the main dipole field can penetrate the solid inner core between infrequent reversals, while [START_REF] Blackett | Magnetic fields of massive rotating bodies[END_REF] proposedanempirical relationship derived from data on the magnetically active bodies in our solar system. He found that their magnetic dipole momentsare proportional to their angular momentum, but thisrule is dimensionally unacceptable because these quantities have unlike units and there is no way to express their ratio in terms of relevant universal constants.

 Magnetic moment:

1/2 5/2 1 Mass Length Time   Angular momentum: 1 2 1
Mass Length Time  However flawed, Blackett's curve-fit may have gotten one key point right. It is entirely possible that ẽxt B  , or equivalently, B TT  , other parameters being equal.

Critique of MAC and IAC balance

A number of authors have attempted to derive scaling rules from assumptions of force or work balance, but many such assumptions lack firm physical foundation. The only rigorous principles of balance are those that derive from three related observations:

 If force terms sum to zero, and if two terms are known to dominate, then those two terms must be roughly equal and opposite. (For example, geostrophic balance assumes equality of pressure gradients and Coriolis forces, but note that such forces do no work.)  If kinetic energy is to remain bounded, then work done on the fluid by buoyancy must over time balance work undone by magnetic and viscous forces.  If the sum of kinetic and magnetic energy is to remain bounded, then the work done by buoyancy must over time balance the energy lost by dissipation as heat.

It is correct but uninformative to equate the work done by buoyancy to the sum of viscous and ohmic dissipation because these diffusive processes operate at vastly shorter length scales, which are not known a-priori.

Magnetic-Archimedean (MA) balance is based on the second observation. It neglects viscosity and equates work done by magnetic (Lorentz) forces to work done by buoyant (Archimedean) forces, which can in turn be related to convective heat transport. The relationship can be made rigorous by adjusting for the ratio of magnetic to total dissipation, which is assumed to be 1 2 () o , supposing rough equipartition of magnetic and kinetic energy in high-order modes and comparable diffusion coefficients:

1 2 [ ( ) ] (1) ( ) PP ag ag r CC o V H o F        V E V B B 
It is a serious (but common) mistake to omit the electric term and to summarize the resulting scaling rule as 22 ~P ag C V B F



. By itself, the magnetic term opposes transverse motions and is inherently dissipative, but in most cases, the electric term offsets or cancels it via induction and/or charge separation. (For example, large-scale flows carry along small-scale magnetic field patterns without feeling an opposing force, and a metal ball can rotate freely in an applied magnetic field parallel to the axis of rotation.)

In other cases of interest, transverse fluid motions distort the magnetic field, but energy is traded elastically, as in Alfvén waves. The concept of magnetic drag is valid only when magnetic diffusion 1/ B    has ample time to smooth the distorted field.

Moreover, diffusion operates at small length scales, whereas the Kolmogorov spectrum shows that buoyant forces do most of their work at large scales.

In convective systems, genuine diffusion is overwhelmed by turbulent transport, so Archimedean-Coriolis (AC) balance asserts that buoyant forces scale with Coriolis forces, but the physical premise of AC balance is suspect because pressure gradients are omitted and Coriolis forces do no work.

-5- The latter prediction combines with MA balance to give 3/10 1/10 BF  , per (Jones 2010).

Lessons from Instability Analysis

Convective instability results from a positive-feedback loop, in which vertical (radial) flows carrywarmer fluid into cooler environments and buoyancy drives said flows. Disregarding horizontal (tangential) components and spatial variation of the flow, in the cavalier spirit of dimensional analysis, we may derive the limiting exponential growth rate by eigenvalue analysis of a linearized system: 0

P d ag dt C rr H HH VV                   



The limiting eigenvalues are

0 P ag C H     
, where 0

d dz H H g     
denotesthe super-adiabaticgradient of enthalpy per unit mass. This relationship is not very useful if the gradient is unknown, but its scaling rule is obvious. Let F denote the heat transported per unit area through a flat convective layer, and Z the distance between the boundaries. Dimensional analysis demands that

3 3 0 ~/ D N F Z   for 1 D N  .
The case where 1 D N  is relevant to the Sun, but data are scarce because in most experiments on heat transport, D N is held constant, while F and ΔTvary. Since it is not possible to convert more than 100% of the heat supplied from below to kinetic energy, we shall use

3 3 0 ~/ FZ   when 1 D N  . We may identify 0 1/ Q T  .
Both Coriolis forces and the main magnetic field act as restraints to radial convection. Stronger temperature gradients are needed to create instability, and in most cases, exponential growth rates are reduced.

-6-Supposing that Coriolis forces convert X Y X V V V    in a flat layer, we could work them into the eigenvalue problem as follows: 00 0 00

P ag d X Z X dt C YY Z H H H VV VV                                When exp( ) i  kx spatial variation is assumed, the relevant eigenvalues are found to be   22 ˆ( ) 1 ( ) ( ) P ag C H          k g k Ωk ,
where hats indicate unit directional vectors.Rapid rotation is seen to suppress the convective instability when k is parallel to Ω, in accord with the Taylor-Proudman theorem. In spherical geometry, rotation preferentially restrains low-M convection modes that vary rapidly with respect to latitude.

There is no satisfactory theory to relate the exponential growth rates inferred from linearized instability analysis to the distribution of energy in convectively driven systems (Barrois 2016a,b). Rigorous dimensional scaling rules à la Kolmogorov could only be stated if all growth rates were rescaled by the same factor.

We cannot conclude that

P ag C H        but 9/5 2/5 2/5 max 0 f      .
By itself, this relationship does not predict V or B, but it can be combined with max Ṽ Z f  and the quasi-rigorous prediction of MA balance:

-7-

2 1 ~P ag C VB FZ  3/ 4 1/ 4 0 3/ 2 1/ 2 0 rms rms B V   

Non-uniform gradients

The gradient is actually non-uniform, so the eigenvalue problem involves a differential equation with boundary condition (0)

( ) 0 V V Z  . 2 2 2 2 ( ) ( ) ( ) ( ) T Z T k V z k H z V z      
One plausible heuristic from the theory of heat transport defines a local turbulent transport coefficient 2 ( ) ẽff zz   , where zrelates to distance from the nearer boundary.

The overall super-adiabatic temperature difference would diverge, were it not for diffusive heat conduction through the boundary layers: .

/2 5/6 0 2 () Z P eff H PH F dz F TF CzC            Actual temperature
If the profile is very steep within boundary layers, but flat elsewhere, then the eigenfunction () Vz can also be relatively flat, except in the boundary layers.

Stratified density

If the density varies with depth, then the equation is altered:

2 2 1 1 2 ( ) ( ) ( ) ( ) T Z Z T k V z k H z V z         
If the density varies very slowly, then () Vz is rescaled by a factor of 1/2   , and the distribution of kinetic energy is essentially unaffected. Near the surface of a star, 2/3 ~~( )

( ) T R r     ,
where  denotes the gravitational potential. The scale height is ill-defined.

Lessons from a high-level dynamo model

An earlier paper by the author [START_REF] Barrois | Bird's-eye view of the solar dynamo[END_REF] proposed a model of the solar dynamo in which three basic processes cooperate to cause oscillations and drive differential rotation. The familiar omega effect bends poloidal field lines () S B into toroidal hairpins () T B , a multi-step regeneration process termed the zeta effect regenerates the dipole field, and their cross-term drives differential rotation () D  . The omega effect depends ondifferential rotation, whereas the zeta effect depends on overall rotation () C  via Coriolis forces on eddies. Linear damping terms and additive noise drivers () Nt  are included to represent the influence of higher-multipole modes omitted from the three-variable model. , which may be correct for the Earth, but oscillatory systems such as the Sun violate the premise of a quasi-constant state.

( ) ( ) ( ) ( ) ( ) ( ) d S C T S S S dt d T D S T T T dt d D S T D D D dt B B g B N t a B B g B N t b B B g N t c                    If
The bird's-eye model predicts that D  will level off provided that 0 ST BB  . A long, wooly, and ultimately inconclusive analysis relates BB to NN , which must in practice be extracted from detailed numerical simulations. The response to a noise process with time dependence exp( ) iwt may be calculated as follows:

1 S S C S T D T T B iw g N B iw g N                               B M N †  BB M NN M
Hypothetically, if NN were simply a multiple of the identity matrix, then one could conclude that There is a further mechanism by which a convectively driven flow mode can pump a weakly damped magnetic field mode, such asa main field component. An -9-analysis of hyperbolic sloshing in (Barrois 2016a) suggests that the amplitude of the pumped magnetic field component should scale with the growth rate of the unstable flow mode.

Comparison with numerical results

The foregoing arguments suggest the following rules for the slow-rotation regime, relevant to the Sun:

 Lorentz number:

  transient fields induced by convective eddies cannot. Recent estimates of conductivity in the outer core (Pozzo 2014) exceed 1.5 MS/m = 1.5 E6 mho/m. Glorified dimensional analysisseeks to identify relevant combinations and to dismiss irrelevant parameters by appealing to physical arguments. It makes short work of three more time scales:  Anything involving the speed of light: / Rc  Anything involving gravity, aside from buoyancy: / Rg  The time for charge separation to cancel electric fields: /

  -Coriolis (IAC) balance asserts that inertial forces scale with buoyant forces,even though inertial forces do no work overall, and also with Coriolis forces. The derivation in[START_REF] Aubert | A systematic experimental study of rapidly rotating spherical convection in water and liquid gallium[END_REF]) can be summarized as follows:

    because a subset of convective modes with k roughly perpendicular to Ω remains unstable. The driven subset comprises a fraction of all modes: token, slowly-varying magnetic fields oppose high-k transverse fluid motions, with 1 () Energy is traded back and forth between the fluid and distorted fields in elastic fashion, setting up Alfvén waves. In spherical geometry, the main toroidal field preferentially restrains high-M convection modes that vary rapidly with longitude, and which are essential to dynamo mechanisms, per Cowling's theorem. The key modes are presumed to have azimuthal wavenumbers of order

  we ignore the noise drivers and eliminate /

  in terms of tensor character. Overall rotation has dipole (L=1) character, whereas differential rotation has octopole (L=3) character.

  But for the rapid-rotation regime relevant to the Earth, the theoretical arguments bog down in unvalidated heuristics, which suggest: poloidal dipole field ext B may scale like rms V , and differential rotation may be inversely related to overall rotation.Christensen&[START_REF] Christensen | Scaling properties of convection-driven dynamos in rotating spherical shells and application to planetary magnetic fields[END_REF] extracted scaling rules from numerical simulations of Earth-like dynamos in the rapid-rotation regime. Little differential rotation was observed, and the poloidal dipole and toroidal quadrupole components were found comparable. rule for the Rossby number agreed with experimental results of Aubert& Brito (2001) and their theoretical rationalization in terms of IAC balance.

Table 1

 1

compares the values of these time scales in the Earth and Sun, but consistent comparisons are difficult because conditions vary with depth and because 1 E-mail: BBarrois@verizon.net -2-many values are uncertain. 2 Sub-surface magnetic fields and fluid flows are largely inaccessible, even though solar flares and helioseismology are revealing.

Table 1 . Comparison of Key Parameters

 1

	Parameter			Formula	Earth	Sun
	Outer radius				outer R		3,480 km	700,000 km
	Inner radius				R		1,220 km	500,000 km
						inner
	Depth of shell	Z R 	R	2,260 km	200,000 km
					outer		inner
	Mass of shell					M		1.15 E23 kg	5.3 E28 kg
	Mean density							11.2 g/cm³	0.058 g/cm³
	Power source					Q		~ 16 E12 watts	3.85 E26 watts
	Polar dipole field					B		~ 5 gauss	~ 1 gauss
						ext
	Differential period				T  		???	10 days
	Rotation period		2/    T	1 day	~ 30 days
	Magnetic energy time scale	B T Z 		/	ext 2 B	17 years	540 years
	Heat transport time scale	Q T		3	2 MZ Q /	1.05 years	0.056 years
	Inductive decay time scale	22 / inner TR    	< 10	4 years	> 10	9 years
	Dissipation number	N	/ DP dr ag C  	< 0.5	6.0

  profiles have been explored experimentally, but the relationship between F and ΔT is shaky. In terms of the Rayleigh number 3 /

	Ra ag T Z  	, and Nusselt number	/  , local Nu FZ C T 		2 ( ) ẽff zz 
		HV							HP
	suggests	1/5 Ñu Ra , whereas global		ẽff Z	2	 suggests	1/2 Ñu Ra . Experiments cited
	by (Faber 1995)have found	0.33 Ñu Ra for	Ra 	8 10	, but	0.28 Ñu Ra for	Ra 	8 10

The magnetic susceptibility, heat capacity, thermal expansion, and transport properties of molten iron mixtures under extreme conditions are all quite uncertain and must be estimated via density functional theory.[START_REF] Anderson | An equation of state for liquid iron and implications for the Earth's core[END_REF][START_REF] Alfè | Iron in the Earth's core: liquid-state thermodynamics and high-pressure melting curve from ab-initio calculations[END_REF][START_REF] Pozzo | Thermal and electrical conductivity of iron and iron-silicon mixtures at Earth's core conditions[END_REF]