Stanis Law

Andrzej Raczyński

Emmanuel Vincent

Shigeki Sagayama

Dynamic Bayesian networks for symbolic polyhonic pitch modeling

The performance of many MIR analysis algorithms, most importantly polyphonic pitch transcription, can be improved by introducing musicological knowledge to the estimation process. We have developed a probabilistically rigorous musicological model that takes into account dependencies between consequent musical notes and consequent chords, as well as the dependencies between chords, notes and the observed note saliences. We investigate its modeling potential by measuring and comparing the cross-entropy with symbolic (MIDI) data.

Introduction

Symbolic pitch modeling, i.e. modeling the prior distribution of note sequences P (N), also known as musicological modeling, is the equivalent of language modeling in speech processing. It has the potential to be used many in Music Information Retrieval (MIR) tasks, like multiple pitch estimation, algorithmic composition, computational musicology, symbolic music analysis, music segmentation, etc., as a part of an integrated statistical model of music [START_REF] Vincent | f(MIR): a roadmap towards versatile MIR[END_REF] .

For example, in the task of polyphonic music transcription, i.e. estimating the pitches, the onset times and the durations of the musical notes present in a recorded audio signal, incorporating a symbolic pitch model would mean a transition from a ML-like estimation N = arg max

N P(S|N) (1)
to estimating the notes in the MAP sense:

N = arg max N P(S|N)P(N), (2)
where P(S|N) is an acoustic model, such as Nonnegative Matrix Factorization 2),3), [START_REF] Vincent | Adaptive harmonic spectral decomposition for multiple pitch estimation[END_REF] , Matching Pursuit, or other model 5),6),7),8),9),?) . While acoustic mod- †1 The University of Tokyo †2 Institut National de Recherche en Informatique et en Automatique eling has been widely studied, symbolic modeling has been given much less attention so far. Some researchers have used basic musicological models in order to overcome the limitations of current state-of-the-art multiple pitch transcription models: Ryynänen and Klapuri 10) proposed a melody transcription method that uses a Hidden Markov Model (HMM) together with a simple musical key model. Their approach however is limited in the sense that it models only monophonic note sequences. Because of that, their approach lacks modeling of the dependencies between concurrent pitches. Raphael and Stoddard [START_REF] Raphael | Harmonic analysis with probabilistic graphical models[END_REF] proposed to use an HMM as a symbolic model for harmonic analysis, i.e. estimating the chord progression behind a sequence of notes. Similar HMMs have also been successfully used for harmonic analysis of audio signals (for a recent paper see e.g. [START_REF] Ueda | HMM-based approach for automatic chord detection using refined acoustic features[END_REF] ). These HMM-based approaches, however, lack absolute pitch modeling and the temporal dependencies are only present between chords.

In this paper we propose a single probabilistic pitch model based on Dynamic Bayesian Networks (DBNs). We model both the dependencies between consequent notes (harmony) and the temporal dependencies between notes and chords. The prior over the note activities P(N) models the temporal dependencies between the hidden variables (similar to those of an HMM) and includes a hidden layer of variables representing chords.

Model

We model the prior distribution of the note sequences P(N) as a DBN with two layers of nodes: the chord (harmony) layer C = (C 1 , C 2 , . . . , C T ) and the note activity layer N = (N 1 , N 2 , . . . , N T ), where T is the number of frames in the analyzed note sequence. We assume that these sequences are first-order Markovian:

P(N) = ∑ C P(C 1 )P(N 1 |C 1 ) • T ∏ t=2 P(N t |N t-1 , C t )P(C t |C t-1 ). ( 3
)
The corresponding network structure is presented in Fig. 1. However, the note activity probability distribution P(N t |N t-1 , C t ) is a highly-dimensional discrete distribution, too complex to train or be used for inference in practice, as it requires 2 88 × 2 88 × 24 (assuming the full piano keyboard) distinct probability values to be fully defined. To deal with this problem, we first factorize the distribution by applying the Bayes rule:

P(N t |N t-1 , C t ) = K ∏ k=1 P(N t,k |N t-1 , C t , N t,1 , N t,2 , . . . , N t,k-1 ), ( 4
)
where k is the analyzed pitch and K = 88 is the size of the analyzed pitch range. This procedure helps to reduce the dimensionality of the note combination variable, but the resulting formula is still difficult to apply in practice due to the highly dimensional conditioning variable set. We therefore split this distribution into three separate ones by using the following approximation: λ3 , where λ = (λ 1 , λ 2 , λ 3 ) are exponential weights, and Z is the normalization factor:

P(N t,k |N t-1 , C t , N t,1 , N t,2 , . . . , N t,k-1 ) (5) ≈Z -1 • P(N t,k |N t-1 ) λ1 • P(N t,k |C t ) λ2 P(N t,k |N t,1 , N t,2 , . . . , N t,k-1 )
Z = 1 ∑ N t,k =0 P(N t,k |N t-1 ) λ1 P(N t,k |C t ) λ2 • (6) • P(N t,k |N t,1 , N t,2 , . . . , N t,k-1 ) λ3 . ( 7
) This kind of log-linear interpolation of musicological models was first proposed by Klakow in the context of language modeling [START_REF] Klakow | Log-linear interpolation of language models[END_REF] . In this paper, each of the separated musicological models is called a submodel and each of them is responsible for modeling a different musicological aspect of the note sequences: voice movement and note duration, harmony and polyphony, respectively.

The exponential weights control the influence of each of the submodels on the resulting note activity distribution: setting them to zero effectively disables the corresponding submodel by making it completely uniform, while using higher values means that the corresponding models will have bigger impact on the results.

Training

The parameters of all submodels were trained by counting occurrences. In all cases a simple smoothing procedure was used: every event was assumed to have occurred at least once.

Chord model

The chord transition probability P(C t |C t-1 ) is easy to model with a multinomial (categorical) probability distribution. This approach is common in MIR tasks that deal with chord progression, e.g. in chord recognition [START_REF] Ueda | HMM-based approach for automatic chord detection using refined acoustic features[END_REF] . It is also common to assume a 24-word chord dictionary, i.e. 12 major and 12 minor chords. We have adopted this approach as well, so the chord transition distribution is described in terms of a 24 × 24 transition matrix.

The left part of Fig. 2 shows the chord transition matrix trained on the entire available dataset. Unfortunately, the obtained transition probabilities are biased, as some keys, and therefore some chord progressions, are sparsely represented in our dataset, while others dominate. However, we can assume that chord transitions have the same distribution in all keys if observed in relation to the tonic, which is reasonable since any song can be transposed to an arbitrary key without any loss in musical correctness. In other words, we assume that the same probability should be given to e.g. the transition from C-major to F-major chord (I→IV transition in C-major key) and the transition from A♭-major to D♭-major (I→IV transition in A♭-major key). In this case, the chord transition probability is a function only of the interval between chord roots and their types. The transition matrix obtained by tying distributions in the above way is presented in Fig. 2.

Furthermore, because key is not considered in our model, we assume a uniform distribution of the initial chord P(C 1 ) = const., which in classical Western music is always the tonic.

Harmony model

Similarly, in order to avoid overfitting, we tie the probabilities of notes having the same musicological function together. This is based on the observation that music can be freely transposed between keys and so all notes should have identical distribution with respect to the chords' root notes. P(N

C C# D D# E F F# G G# A A# B c c# d d# e f f# g g# a a# b C C# D D# E F F# G G# A A# B c c# d d# e f f# g g# a a# b M->M m->M m->M m->m Chord transition matrix Previous chord Current chord C C# D D# E F F# G G# A A# B c c# d d# e f f# g g# a a# b C C# D D# E F F# G G# A A# B c c# d d# e f f# g g# a a# b M->M m->M m->M m->m
Fig. mode is the mode operator, i.e. major or minor. The corresponding probability distribution is presented in Fig. 3.

Voice model

The voice model P(N t,k |N t-1 ) also suffers from the high dimensionality of the conditioning variable set. To reduce the dimensionality, we have proposed three different, simplified models that assume that only some pitches are relevant in modeling the voice movement.

Reduced voice model

In this model we assume dependence only on pitches from a small range (-R, . . . , 0, 1, . . . , R), relative to the current pitch k:

P(N t,k |N t-1 ) ≈ P(N t,k |N t-1,k-R , N t-1,k-R+1 , . . . , N t-1,k , . . . , N t-1,k+R ), (9)
This is based on an assumption that voice movement is limited to small jumps, typically within a single octave, i.e. R = 12.

Duration-only model

In this model individual note activities are assumed to be dependent only on the previous state of the same pitch:

P(N t,k |N t-1 ) ≈ P(N t,k |N t-1,k ), (10)
i.e. by a pitch-dependent conditional Bernoulli model.

Voice movement model

In this model we use the following approximation:

P(N t,k |N t-1 ) ≈ P(N t,k |M t,k ), (11)
where M t,k is the distance between the pitch k and the closest active pitch in the previous time frame. If the pitch k was active in the previous time frame, this model acts as a duration model, otherwise it is a simple voice movement model.

Comparison

The models are compared by the cross-entropies of the observed note data N (see section 5), calculated on the testing data set (see section 4). We have compared the duration-only model (D), the reduced voice (R) and the durationonly model models combined by log-linear interpolation (R+D), and two versions of the voice movement model: independent voice movement model (V) and an independent version identically distributed for all pitches (iid. V). For a baseline we have juxtaposed these models with an unconditional Bernoulli model (B). The results are presented in Fig. 4. ). The dashed line marks the marginal note activity probability P(N t,k ).

0

Polyphony model

The probability distribution P(N t,k |N t,1 , N t,2 , . . . , N t,k-1 ) is more difficult to model due to high dimensionality of the conditioning variable set. However, since global note distributions and pitch are already modeled by the duration and the harmony model, respectively, the polyphony model can be simplified to model only the number of notes active simultaneously:

P(N t,k |N t,1 , N t,2 , . . . , N t,k-1 ) ≈ P(N t,k |L t,k ), ( 12
) where L t,k = ∑ k-1 m=1 N t,m
. The resulting distribution is plotted in Fig. 5.

Exponential weights

The interpolation weights λ from Eq. 5 are optimized by maximizing their log-likelihood: λ = arg max

λ log P(N|λ), (13)
which is calculated for the validation dataset (see section 4). Optimization is performed using the Broyden-Fletcher-Goldfarb-Shanno (BFGS) method (a quasi-Newton optimization), built in the GNU R environment as the optim() function [START_REF]R: A Language and Environment for Statistical Computing[END_REF] . The resulting values are listed in values were all equal to λ p = 1, p ∈ {1, 2, 3}.

Data

Two datasets were used in the experiments: the widely used RWC database [START_REF] Goto | RWC music database: Music genre database and musical instrument sound database[END_REF] and the Mutopia Project dataset [START_REF]The Mutopia Project. Free classical and contemporary sheet music[END_REF] . The classical pieces of the RWC database were annotated with detailed harmony labels that include: keys and modulations, and chords with their roots, inversions, types and various modifications [START_REF] Kaneko | Functional harmony annotation database for statistical music analysis[END_REF] . This data uses abstract, tempo-independent musical time (measures and beats), and served as the chord ground-truth for training the harmony and chord models.

The Mutopia dataset consisted of 1468 files divided into 3 subsets: for training (1268 files), validation (100 files) and testing (100 files). The training set was used to train all remaining submodels, while the submodel weights from Eq. 5 were trained on the validation set. Experiments were performed on the testing data. All symbolic data was quantized and 1 frame corresponded to 1/6 th of a beat, so it was rather sheet music encoded in MIDI format than performance data, as is the case with the original RWC data.

Symbolic evaluation

The models are compared by calculating the cross-entropy of the observed note data given the musicological prior. We compute the marginal crossentropies, i.e. cross-entropies with the chord sequence marginalized out:

H(N) = - 1 T log 2 P(N|λ) (14) = - 1 T log 2 ∑ CO P(N|C O , λ)P(C O |λ).
The conditional cross-entropy H(N|C) is calculated by summing over all possible chord sequences. This is calculated using the frontier algorithm [START_REF] Murphy | Dynamic Bayesian networks: representation, inference and learning[END_REF] , which is the DBN equivalent of the forward-backward algorithm. Three reference models were used for comparison: the uniform model P(N t,k = 1) = 1/2, the The results are presented in Fig. 7. They how that modeling harmonic relations between pitches (HC) leads to better modeling performance than modeling prior pitch distributions with a Bernoulli model (i.). However, modeling horizontal relations between pitches (V) seems more to be important than the vertical harmony models. Combining all proposed model (HCVP) results in the lowest note entropy.

Conclusion and future work

We have presented a probabilistic symbolic pitch model that is a log-linear interpolation of a number of simpler models representing complementary properties of a note sequence. The proposed models have been evaluated by calculating the marginal cross-entropy H(N) of the testing data set given the model. The entropy of the observed notes is as low as 6.44 bits per time frame, compared to the reference unconditional Bernoulli model (18.03 bits) and the uniform model (88 bits). In future work we will focus on combining the proposed models with other MIR models from 1) , e.g. an acoustic model to perform multiple pitch estimation, so as to observe the practical implications of incorporating a musicological prior to the estimation process.

Acknowledgments

This work is supported by INRIA under the Associate Team Program VER-SAMUS (http://versamus.inria.fr/).

Fig. 1

 1 Fig. 1 Proposed structure of the Dynamic Bayesian Network for polyphonic pitch modeling.

Fig. 3

 3 Fig. 3 Pitch probability distribution for major (M) and minor (m) chords as a function of the interval from the chord's root note.

Fig. 4

 4 Fig. 4 Cross-entropy obtained for different voice models and an unconditional Bernoulli model (B) for reference.

Fig. 6

 6 Fig. 6 Parameters p k of the independent note activity model. Black and white bars correspond to black and white piano keys, respectively.

Fig. 7

 7 Fig. 7 Average marginal cross-entropies for the tested models. iid. = iid. Bernoulli, i. = independent Bernoulli, H = Harmony model, C = chord model, P = Polyphony model, V = Voice model.

  t,k |C t ) = P( inter{k; root{C t }}|mode{C t }),(8) where inter is the musical interval operator, root is the root note operator and

		Vol.2011-MUS-91 No.8
		2011/7/28
	Chord transition matrix
	Current chord
	Previous chord
	2	c ⃝ 2011 Information Processing Society of Japan

Table 1 .

 1 The initial coefficient

	Vol.2011-MUS-91 No.8
	2011/7/28

c ⃝ 2011 Information Processing Society of Japan 4

Table 1

 1 Trained values of exponential coefficients for submodels from Eq. 5.


Ahmed Bouajjani

Michael Emmi

Bounded Phase Analysis of Message-Passing Programs *

We describe a novel technique for bounded analysis of asynchronous message-passing programs with ordered message queues. Our bounding parameter does not limit the number of pending messages, nor the number of "contexts-switches" between processes. Instead, we limit the number of process communication cycles, in which an unbounded number of messages are sent to an unbounded number of processes across an unbounded number of contexts. We show that remarkably, despite the potential for such vast exploration, our bounding scheme gives rise to a simple and efficient program analysis by reduction to sequential programs. As our reduction avoids explicitly representing message queues, our analysis scales irrespectively of queue content and variation.

Introduction

Software is becoming increasingly concurrent: reactivity (e.g., in user interfaces, web servers), parallelization (e.g., in scientific computations), and decentralization (e.g., in web applications) necessitate asynchronous computation. Although sharedmemory implementations are often possible, the burden of preventing unwanted thread interleavings without crippling performance is onerous. Many have instead adopted asynchronous programming models in which processes communicate by posting messages/tasks to others' message/task queues-Miller et al. [START_REF] Miller | Concurrency among strangers[END_REF] discuss why such models provide good programming abstractions. Single-process systems such as the JavaScript page-loading engine of modern web browsers [START_REF]HTML5: A vocabulary and associated APIs for HTML and XHTML[END_REF], and the highly-scalable Node.js asynchronous web server [START_REF] Dahl | Node.js: Evented I/O for V8 JavaScript[END_REF], execute a series of short-lived tasks one-by-one, each task potentially queueing additional tasks to be executed later. This programming style ensures that the overall system responds quickly to incoming events (e.g., user input, connection requests). In the multi-process setting, languages such as Erlang and Scala have adopted messagepassing as a fundamental construct with which highly-scalable and highly-reliable distributed systems are built.

Despite the increasing popularity of such programming models, little is known about precise algorithmic reasoning. This is perhaps not without good reason: decision problems such as state-reachability for programs communicating with unbounded reliable queues are undecidable [START_REF] Brand | On communicating finite-state machines[END_REF], even when there is only a single finite-state process (posting messages to itself). Furthermore, the known decidable under-approximations (e.g., bounding the size of queues) represent queues explicitly, are thus doomed to combinatorial explosion as the size and variability of queue content increases. Some have proposed analyses which abstract message arrival order [START_REF] Sen | Model checking multithreaded programs with asynchronous atomic methods[END_REF][START_REF] Jhala | Interprocedural analysis of asynchronous programs[END_REF][START_REF] Ganty | Algorithmic verification of asynchronous programs[END_REF], or assume messages can be arbitrarily lost [START_REF] Abdulla | Verifying programs with unreliable channels[END_REF][START_REF] Abdulla | On-the-fly analysis of systems with unbounded, lossy fifo channels[END_REF]. Such analyses do not suffice when correctness arguments rely on reliable messaging-several systems specifically do ensure the ordered delivery of messages, including Scala, and recent web-browser specifications [START_REF]HTML5: A vocabulary and associated APIs for HTML and XHTML[END_REF]. Others have proposed analyses which compute finite symbolic representations of queue contents [START_REF] Boigelot | Symbolic verification of communication protocols with infinite state spaces using QDDs[END_REF][START_REF] Bouajjani | Verification of parametric concurrent systems with prioritised FIFO resource management[END_REF]. Known bounded analyses which model queues precisely either bound the maximum capacity of message-queues, ignoring executions which exceed the bound, or bound the total number of process "contexts" [START_REF] Qadeer | Context-bounded model checking of concurrent software[END_REF][START_REF] La Torre | Context-bounded analysis of concurrent queue systems[END_REF], where each context involves a single process sending and receiving messages. For each of these bounding schemes there are trivial systems which cannot be adequately explored, e.g., by sending more messages than the allowed queue-capacity, having more processes than contexts, or by alternating message-sends to two processes-we discuss such examples in Section 3. All of the above techniques represent queues explicitly, though perhaps symbolically, and face combinatorial explosion as queue content and variation increase.

In this work we propose a novel technique for bounded analysis of asynchronous message-passing programs with reliable, ordered message queues. Our bounding parameter, introduced in Section 3, is not sensitive to the capacity nor content of message queues, nor the number of process contexts. Instead, we bound the number of process communication cycles by labeling each message with a monotonically-increasing phase number. Each time a message chain visits the same process, the phase number must increase. For a given parameter k, we only explore behaviors of up to k phases-though k phases can go a long way. In the leader election distributed protocol [START_REF] Svensson | A new leader election implementation[END_REF] for example, each election round occurs in 2 phases: in the first phase each process sends capture messages to the others; in the second phase some processes receive accept messages, and those that find themselves majority-winners broadcast elected messages. In these two phases an unbounded number of messages are sent to an unbounded number of processes across an unbounded number of process contexts! We demonstrate the true strength of phase-bounding by showing in Sections 4 and 5 that the bounded phase executions of a message-passing program can be concisely encoded as a non-deterministic sequential program, in which messagequeues are not explicitly represented. Our so-called "sequentialization" sheds hope for scalable analyses of message-passing programs. In a small set of simple experiments (Section 4), we demonstrate that our phase-bounded encoding scales far beyond known explicit-queue encodings as queue-content increases, and even remains competitive as queue-content is fixed while the number of phases grows. By reducing to sequential programs, we leverage highly-developed sequential program analysis tools for message-passing programs.

Asynchronous Message-Passing Programs

We consider a simple multi-processor programming model in which each processor is equipped with a procedure stack and a queue of pending tasks. Initially all processors are idle. When an idle processor's queue is non-empty, the oldest task in its queue is removed and executed to completion. Each task executes essentially a recursive sequential program, which besides accessing its own processor's global storage, can post tasks to the queues of any processor, including its own. When a task does complete, its processor again becomes idle, chooses the next pending task to execute to completion, and so on. The distinction between queues containing messages and queues containing tasks is mostly aesthetic, but in our task-based treatment queues are only read by idle processors; reading additional messages during a task's execution is prohibited. While in principle many message-passing systems, e.g., in Erlang and Scala, allow reading additional messages at any program point, we have observed that common practice is to read messages only upon completing a task [START_REF] Trottier-Hebert | Learn you some Erlang for great good![END_REF].

Though similar to Sen and Viswanathan [START_REF] Sen | Model checking multithreaded programs with asynchronous atomic methods[END_REF]'s model of asynchronous programs, the model we consider has two important distinctions. First, tasks execute across potentially several processors, rather than only one, each processor having its own global state and pending tasks. Second, the tasks of each processor are executed in exactly the order they are posted. For the case of single-processor programs, Sen and Viswanathan [START_REF] Sen | Model checking multithreaded programs with asynchronous atomic methods[END_REF]'s model can be seen as an abstraction of the model we consider, since there the task chosen to execute next when a processor is idle is chosen non-deterministically among all pending tasks.

Program Syntax

Let Procs be a set of procedure names, Vals a set of values, Exprs a set of expressions, Pids a set of processor identifiers, and let T be a type. Figure 1 gives the grammar of asynchronous message-passing programs. We intentionally leave the syntax of expressions e unspecified, though we do insist Vals contains true and false, and Exprs contains Vals and the (nullary) choice operator .

Each program P declares a single global variable g and a procedure sequence, each p ∈ Procs having a single parameter l and top-level statement denoted s p ; as statements are built inductively by composition with control-flow statements, s p describes the entire body of p. The set of program statements s is denoted Stmts. Intuitively, a post ρ p e statement is an asynchronous call to a procedure p with argument e to be executed on the processor identified by ρ; a self-post to one's own processor is made by setting ρ to . A program in which all post statements are self-posts is called a single-processor program, and a program without post statements is called a sequential program.

The programming language we consider is simple, yet very expressive, since the syntax of types and expressions is left free, and we lose no generality by considering only single global and local variables. Appendix A lists several syntactic extensions which we use in the source-to-source translations of the subsequent sections, and which easily reduce to the syntax of our grammar. The grammar of asynchronous message-passing programs P . Here T is an unspecified type, and e, p, and ρ range, resp., over expressions, procedure names, and processor identifiers.

Dispatch g, ε, f • q S -→ g, f, q Complete f = , return e; s g, f, q S -→ g, ε, q Self-Post s1 = post p e; s2
2 ∈ e(g, 1) f = 2, sp g, 1, s1 w, q S -→ g, 1, s2 w, q • f Fig. 2. The single-processor transition rules → S ; see Appendix B for the standard sequential statements.

Single-Processor Semantics

A (procedure) frame f = , s is a current valuation ∈ Vals to the procedurelocal variable l, along with a statement s ∈ Stmts to be executed. (Here s describes the entire body of a procedure p that remains to be executed, and is initially set to p's top-level statement s p ; we refer to initial procedure frames t = , s p as tasks, to distinguish the frames that populate processor queues.)

The set of all frames is denoted Frames.

A processor configuration κ = g, w, q is a current valuation g ∈ Vals to the processor-global variable g, along with a procedure-frame stack w ∈ Frames * and a pending-tasks queue q ∈ Frames * . A processor is idle when w = ε. The set of all processor configurations is denoted Pconfigs. A processor configuration map ξ : Pids → Pconfigs maps each processor ρ ∈ Pids to a processor configuration ξ(ρ). We write ξ (ρ → κ) to denote the configuration ξ updated with the mapping (ρ → κ), i.e., the configuration ξ such that ξ (ρ) = κ, and ξ (ρ ) = ξ(ρ ) for all ρ ∈ Pids \ {ρ}.

For expressions without program variables, we assume the existence of an evaluation function • e : Exprs → ℘(Vals) such that e = Vals. For convenience we define e(g, ) def = e[g/g, /l] e to evaluate the expression e in a global valuation g by substituting the current values for variables g and l. As these are the only program variables, the substituted expression e[g/g, /l] has no free variables.

Figure 2 defines the transition relation → S for the asynchronous behavior of each processor; the standard transitions for the sequential statements are listed in Appendix B. The Self-Post rule creates a new frame to execute the given procedure, and places the new frame in the current processor's pending-tasks queue. The Complete rule returns from the final frame of a task, rendering the processor idle, and the Dispatch rule schedules the least-recently posted task on a idle processor.

Switch ρ2 ∈ enabled(m, ξ) ρ1, ξ, m -→ M ρ2, ξ, m Step ξ1(ρ) S -→ κ ξ2 = ξ1 (ρ → κ) ρ ∈ enabled(m1, ξ1) m2 = step(m1, ξ1, ξ2) ρ, ξ1, m1 -→ M ρ, ξ2, m2
Post ξ1(ρ1) = g1, 1, post ρ2 p e; s w1, q1 ξ1(ρ2) = g2, w2, q2 ρ1 = ρ2

2 ∈ e(g1,

1) f = 2, sp ρ1 ∈ enabled(m1, ξ1) m2 = step(m1, ξ1, ξ3) ξ2 = ξ1 (ρ1 → g1, 1, s w1, q1 ) ξ3 = ξ2 (ρ2 → g2, w2, q2 • f ) ρ1, ξ1, m1 -→ M ρ1, ξ3, m2

Multi-Processor Semantics

In reality the processors of multi-processor systems execute independently in parallel. However, as long as they either do not share memory, or access a sequentially consistent shared memory, it is equivalent, w.r.t. the observations of any single processor, to consider an interleaving semantics: at any moment only one processor executes. In order to later restrict processor interleaving, we make explicit the scheduler which arbitrates the possible interleavings. there is some ρ ∈ Pids such that ξ(ρ) is either non-idle or has pending tasks, then there exists ρ ∈ Pids such that ρ ∈ enabled(m, ξ) and ξ(ρ ) is either non-idle or has pending tasks. A configuration c = ρ, ξ, m is a currently executing processor ρ ∈ Pids, along with a processor configuration map ξ, and a scheduler object m. Figure 3 defines the multi-processor transition relation → M , parameterized by a scheduler M . The Switch rule non-deterministically schedules any enabled processor, while the Step rule executes one single-processor program step on the currently scheduled processor, and updates the scheduler object. Finally, the Post rule creates a new frame to execute the given procedure, and places the the new frame on the target processor's pending-tasks queue.

Until further notice, we assume M is a completely non-deterministic scheduler; i.e., all processors are always enabled. In Section 5 we discuss alternatives.

An M -execution of a program P (from c 0 to c j ) is a configuration sequence c 0 c 1 . . . c j such that c i → M c i+1 for 0 ≤ i < j. An initial condition ι = ρ 0 , g 0 , 0 , p 0 is a processor identifier ρ 0 , along with a global-variable valuation g 0 ∈ Vals, a local-variable valuation 0 ∈ Vals, and a procedure p 0 ∈ Procs. A configuration c = ρ 0 , ξ, empty of a program P is ρ 0 , g 0 , 0 , p 0 -initial when ξ(ρ 0 ) = g 0 , ε, 0 , s p0 and ξ(ρ) = g 0 , ε, ε for all ρ = ρ 0 . A configuration ρ, ξ, m is g f -final when ξ(ρ ) = g f , w, q for some ρ ∈ Pids, and w, q ∈ Frames * . We say a global valuation g is M -reachable in P from ι when there exists an M -execution of P from some c 0 to some c such that c 0 is ι-initial and c is g-final 1 .

Definition 1. The state-reachability problem is to determine for an initial condition ι, valuation g, and program P , whether g is reachable in P from ι.

Phase-Bounded Execution

Because processors execute tasks precisely in the order which they are posted to their unbounded task-queues, our state-reachability problem is undecidable, even with only a single processor accessing finite-state data [START_REF] Brand | On communicating finite-state machines[END_REF]. Since it is not algorithmically possible to consider every execution precisely, in what follows we present an incremental under-approximation. For a given bounding parameter k, we consider a subset of execution (prefixes) precisely; as k increases, the set of considered executions increases, and in the limit as k approaches infinity, every execution of any program is considered-though for many programs, every execution is considered with a finite value of k.

In a given execution, a task-chain t 1 t 2 . . . t i from t 1 to t i is a sequence of tasks2 such that the execution of each t j posts t j+1 , for 0 < j < i, and we say that t 1 is an ancestor of t i . We characterize execution prefixes by labeling each task t posted in an execution with a phase number ϕ(t) ∈ N:

ϕ(t) =          0 if t is initially pending. ϕ(t )
if t is posted to processor ρ by t , and t has no phase-ϕ(t ) ancestor on processor ρ.

ϕ(t ) + 1 if t is posted by t , otherwise.
For instance, considering Figure 5a, supposing all on a single processor, an initial task A 1 posts A 2 , A 3 , and A 4 , then A 2 posts A 5 and A 6 , and then A 3 posts A 7 , which in turn posts A 8 and A 9 . Task A 1 has phase 0. Since each post is made to the same processor, the phase number is incremented for each posted task. Thus the phase 1 tasks are {A 2 , A 3 , A 4 }, the phase 2 tasks are {A 5 , A 6 , A 7 }, and the phase 3 tasks are {A 8 , A 9 }. Notice that tasks of a given phase only execute after all tasks of the previous phase have completed, i.e., execution order is in phase order; only executing tasks up to a given phase does correspond to a valid execution prefix.

Definition 2. An execution is k-phase when ϕ(t) < k for each executed task t.

The execution in Figure 5a is a 4-phase execution, since all tasks have phase less than 4. Despite there being an arbitrary number 4n + 1 of posted tasks, the execution in Figure 5b is 1-phase, since there are no task-chains between same-processor tasks. Contrarily, the execution in Figure 5c requires n phases to execute all 2n tasks, since every other occurrence of an A i task creates a task-chain between A-tasks.

Note that bounding the number of execution phases does not necessarily bound the total number of tasks executed, nor the maximum size of task queues, nor the amount of switching between processors. Instead, a bound k restricts the maximum length of task chains to k • |Pids|. In fact, phase-bounding is incomparable to bounding the maximum size of task queues. On the one hand, every execution of a program in which one root task posts an arbitrary, unbounded number of tasks to other processors (e.g., in Figure 5b) are explored with 1 phase, though no bound on the size of queues will capture all executions. On the other hand, all executions with a single arbitrarily-long chain of tasks (e.g., in Figure 5c) are explored with size 1 task queues, though no limited number of phases captures all executions. In the limit as the bounding parameter increases, either scheme does capture all executions.

Theorem 1 (Completeness). For every execution h of a program P , there exists k ∈ N such that h is a k-phase execution.

Phase-Bounding for Single-Processor Programs

Characterizing executions by their phase-bound reveals a simple and efficient technique for bounded exploration. This seems remarkable, given that phase-bounding explores executions in which arbitrarily many tasks execute, making the task queue arbitrarily large. The first key ingredient is that once the number of phases is bounded, each phase can be executed in isolation. For instance, consider again the execution of Figure 5a. In phase 1, the tasks A 2 , A 3 , and A 4 pick up execution from the global valuation g 1 which A 1 left off at, and leave behind a global valuation g 2 for the phase 2 tasks. In fact, given the sequence of tasks in each phase, the only other "communication" between phases is a single passed global valuation; executing that sequence of tasks on that global valuation is a faithful simulation of that phase.

The second key ingredient is that the ordered sequence of tasks executed in a given phase is exactly the ordered sequence of tasks posted in the previous phase. This is obvious, since tasks are executed in the order they are posted. However, combined with the first ingredient we have quite a powerful recipe. Supposing the global state g i at the beginning of each phase i is known initially, we can simulate a k-phase execution by executing each task posted to phase i as soon as it is posted, with an independent virtual copy of the global state, initially set to g i . That is, our simulation will store a vector of k global valuations, one for each phase. Initially, the i th global valuation is set to the state g i in which phase i begins; tasks of phase i then read from and write to the i th global valuation. It then only remains to ensure that the global valuations g i used at the beginning of each phase 0 < i < k match the valuations reached at the end of phase i -1.

This simulation is easily encoded into a non-deterministic sequential program with k copies of global storage. The program begins by non-deterministically setting each copy to an arbitrary value. Each task maintains their current phase number i, and accesses the i th copy of global storage. Each posted task is simply called instead of posted, its phase number set to one greater than its parentposts to tasks with phase number k are ignored. At the end of execution, the program ensures that the i th global valuation matches the initially-used valuation for phase i + 1, for 0 ≤ i < k -1. When this condition holds, any global valuation observed along the execution is reachable within k phases in the original program. Figure 4 lists a code-to-code translation which implements this simulation. Theorem 2. A global-valuation g is reachable in a k-phase execution of a singleprocessor program P if and only if g is reachable in ((P )) k -the k-phase sequential translation of P .

When the underlying sequential program model has a decidable state-reachability problem, Theorem 2 gives a decision procedure for the phase-bounded statereachability problem, by applying the decision procedure for the underlying model to the translated program. This allows us for instance to derive a decidability result for programs with finite data domains.

Corollary 1. The k-phase state-reachability problem is decidable for singleprocessor programs with finite data domains.

More generally, given any underlying sequential program model, our translation makes applicable any analysis tool for said model to message-passing programs, since the values of the additional variables are either from the finite domain {0, . . . , k -1}, or in the domain of the original program variables.

Note that our simulation of a k-phase execution does not explicitly store the unbounded task queue. Instead of storing a multitude of possible unbounded task sequences, our simulation stores exactly k global state valuations. Accordingly, our simulation is not doomed to the unavoidable combinatorial explosion encountered by storing (even bounded-size) task queues explicitly. To demonstrate the capability of our advantage, we measure the time to verify two fabricated yet illustrative examples (listed in full in Appendix C, comparing our bounded-phase encoding with a bounded task-queue encoding. In the bounded task-queue encoding, we represent the task-queue explicitly by an array of integers, which stores the identifiers of posted procedures 3 . When control of the initial task completes, the program enters a loop which takes a procedure identifier from the head of the queue, and calls the associated procedure. When the queue reaches a given bound, any further posted tasks are ignored.

The first program P 1 (i), parameterized by i ∈ N, has a single Boolean global variable b, i procedures named p 1 , . . . , p i , which assert b to be false and set b to true, and i procedures named q 1 , . . . , q i which set b to false. Initially, P 1 (i) sets b to false, and enters a loop in which each iteration posts some p j followed by some q j . Since a q j task must be executed between each p j task, each of the assertions are guaranteed to hold. Figure 6a compares the time required to verify P 1 (i) (using the Boogie verification engine [START_REF] Barnett | Boogie: An intermediate verification language[END_REF]) for various values of i, and various bounds n on loop unrolling. Note that although every execution of P 1 (i) has only 2 phases, to explore all n loop iterations in any given execution, the size of queues must be at least 2n, since two tasks are posted per iteration. Even for this very simple program, representing (even bounded) task-queues explicitly does not scale, since the number of possible task-queues grows astronomically as the size of task-queues grow. This ultimately prohibits the bounded tasks-queue encodings from exploring executions in which more than a mere few simple tasks execute. On the contrary, our bounded-phase simulation easily explores every execution up to the loop-unrolling bound in a few seconds.

To be fair, our second program P 2 is biased to support the bounded task-queue encoding. Following the example of Figure 5c, P 2 again has a single Boolean global variable b, and two procedures: p 1 asserts b to be false, sets b to true, and posts p 2 , while p 2 sets b to false and posts p 1 . Initially, the program P 2 sets b to false and posts a single p 1 task. Again here, since a p 2 task must execute between each p 1 task, each of the assertions are guaranteed to hold. Figure 6b compares the time required to verify P 2 for various bounds n on the number of tasks explored 4 . Note that although every execution of P 2 uses only size 1 task-queues, to explore all n tasks in any given execution, the number of phases 3 For simplicity our examples do not pass arguments to tasks; in general, one should also store in the task-queue array the values of arguments passed to each posted procedure. 4 The number n of explored tasks is controlled by limiting the number of loop unrollings in the bounded task-queue encoding, and limiting the recursion depth, and phase-bound, in the bounded-phase encoding. must be at least n, since each task must execute in its own phase. Although verification time for the bounded-phase encoding does increase with n faster than the bounded task-queue encoding-as expected-due to additional copies of the global valuation, and more deeply in-lined procedures, the verification time remains manageable. In particular, the time does not explode uncontrollably: even 50 tasks are explored in under 20s.

Phase-Bounding for Multi-Processor Programs

Though state-reachability under a phase bound is immediately and succinctly reducible to sequential program analysis for single-processor programs, the multiprocessor case is more complicated. The added complexity arises due to the many orders in which tasks on separate processors can contribute to others' task-queues. As a simple example, consider the possible bounded-phase executions of Figure 5b with four processors, A, B, C, and D. Though B's tasks B 1 , . . . , B n must be executed in order, and C's tasks C 1 , . . . , C n must also be executed in order, the order of D's tasks are not pre-determined: the arrival order of D's tasks depends on how B's and C's tasks interleave. Suppose for instance B 1 executes to completion before C 1 , which executes to completion before B 2 , and so on. In this case D's tasks arrive to D's queue, and ultimately execute, in the index order D 1 , D 2 , . . . as depicted. However, there exist executions for every possible order of D's tasks respecting D 1 < D 3 < . . . and D 2 < D 4 < . . . (where < denotes an ordering constraint)-many possible orders indeed! In fact, due to the capability of such unbounded interleaving, the problem of state-reachability under a phase-bound is undecidable for multi-processor programs, even for programs with finite data domains.

Theorem 3. The k-phase bounded state-reachability problem is undecidable for multi-processor programs with finite data domains.

Note that Theorem 3 holds independently of whether memory is shared between processors: the fact that a task-queue can store any possible (unbounded) shuffling of tasks posted by two processors lends the power to simulate Post's correspondence problem [START_REF] Post | A variant of a recursively unsolvable problem[END_REF].

Theorem 3 insists that phase-bounding alone will not lead to the elegant encoding to sequential programs which was possible for single-processor programs. If that were possible, then the translation from a finite-data program would lead to a finite-data sequential program, and thus a decidable state-reachability problem. Since a precise algorithmic solution to bounded-phase state-reachability is impossible for multi-processor programs, we resort to a further incremental yet orthogonal under-approximation, which limits the number of considered processor interleavings. The following development is based on delay-bounded scheduling [START_REF] Emmi | Delay-bounded scheduling[END_REF].

We define a delaying scheduler M = D, empty, enabled, step, delay , as a scheduler D, empty, enabled, step , along with a function delay : (D × Pids × (Pids → Pconfigs)) → D. Furthermore, we extend the transition relation of Figure 3 with a postponing rule of Figure 7 which we henceforth refer to as a delay (operation), saying that processor ρ is delayed. Note that a delay operation may or may not change the set of enabled processors in any given step, depending on the scheduler. A delaying scheduler is delay-accessible when for every configuration c 1 and non-idle or task-pending processor ρ, there exists a sequence c 1 → M . . . → M c j of Delay-steps such that ρ is enabled in c j . Given executions h 1 and h 2 of (delaying) schedulers M 1 and M 2 resp., we write h 1 ∼ h 2 when h 1 and h 2 are identical after projecting away delay operations. Definition 3. An execution with at most k delay operators is called k-delay.

Consider again the possible executions of Figure 5b, but suppose we fix a deterministic scheduler M which without delaying would execute D's tasks in index order: D 1 , D 2 , . . .; furthermore suppose that delaying a processor ρ in phase i causes M to execute the remaining phase i tasks of ρ in phase i + 1, while keeping the tasks of other processors in their current phase. Without using any delays, the execution of Figure 5b is unique, since M is deterministic. However, as Figure 8 illustrates, using a single delay, it is possible to also derive the order D 1 , D 3 , . . . , D 2n-1 , D 2 , D 4 , . . . , D 2n (among others): simply delay processor C once before C 1 posts D 2 . Since this forces the D 2i tasks posted by each C i to occur in the second phase, it follows they must all happen after the D 2i-1 tasks posted by each B i .

Theorem 4 (Completeness). Let M be any delay-accessible scheduler. For every execution h of a program P , there exists an M -execution h and k ∈ N such that h is a k-delay execution and h ∼ h.

Note that Theorem 4 holds for any delay-accessible scheduler M -even deterministic schedulers. As it turns out there is one particular scheduler M bfs for which we know a convenient sequential encoding, and this scheduler is described in Appendix D. For the moment, the important points to note are that M bfs is deterministic, non-blocking, and delay-accessible. Essentially, determinism allows us to encode the scheduler succinctly in a sequential program; the nonblocking property ensures this scheduler does explore some execution, rather than needlessly ceasing to continue; delay-accessibility combined with Theorem 4 ensure the scheduler is complete in the limit. Figure 9 lists a code-to-code translation which encodes bounded-phase and bounded-delay exploration of a given program according to the M bfs scheduler as a sequential program.

Our translation closely follows the single-processor translation of Section 4, the key differences being:

the phase of a posted task is not necessarily incremented, since posted tasks may not have same-processor ancestors in the current phase, and at any point, the currently executing task may increment a delay counter, causing all following tasks on the same processor to shift forward one additional phase.

As the global values reached by each processor at the end of each phase i -1 must be ensured to match the initial values of phase i, for 0 < i < k + d, so must the values for the shift counter: an execution is only valid when for each processor ρ ∈ Pids and each phase 0

< i < k, shift[ρ][i -1] matches the initial value of shift[ρ][i].
Theorem 5. A global valuation g is reachable in a k-phase d-delay M bfs -execution of a multi-processor program P if and only if g is reachable in

((P )) bfs k,d .
As is the case for our single-processor translation, our simulation does not explicitly store the unbounded tasks queue, and is not doomed to combinatorial explosion faced by storing tasks-queues explicitly.

Related Work

Our work follows the line of research on compositional reductions from concurrent to sequential programs. The initial so-called "sequentialization" [START_REF] Qadeer | KISS: Keep it simple and sequential[END_REF] explored multi-threaded programs up to one context-switch between threads, and was later expanded to handle a parameterized amount of context-switches between a statically-determined set of threads executing in round-robin order [START_REF] Qadeer | Context-bounded model checking of concurrent software[END_REF][START_REF] Lal | Reducing concurrent analysis under a context bound to sequential analysis[END_REF]. La Torre et al. [START_REF] La Torre | Model-checking parameterized concurrent programs using linear interfaces[END_REF] later extended the approach to handle programs parameterized by an unbounded number of statically-determined threads, and shortly after, Emmi et al. [START_REF] Emmi | Delay-bounded scheduling[END_REF] further extended these results to handle an unbounded amount of dynamically-created tasks, which besides applying to multi-threaded programs, naturally handles asynchronous event-driven programs [START_REF] Sen | Model checking multithreaded programs with asynchronous atomic methods[END_REF]. Bouajjani et al. [START_REF] Bouajjani | On sequentializing concurrent programs[END_REF] pushed these results even further to a sequentialization which attempts to explore as many behaviors as possible within a given analysis budget. Each of these sequentializations necessarily do provide a bounding parameter which limits the amount of interleaving between threads or tasks, but none are capable of precisely exploring tasks in creation order, which is abstracted away from their program models [START_REF] Sen | Model checking multithreaded programs with asynchronous atomic methods[END_REF]. Kidd et al. [START_REF] Kidd | One stack to run them all: Reducing concurrent analysis to sequential analysis under priority scheduling[END_REF]'s sequentialization is sensitive to task priorities, their reduction assumes a finite number of statically-determined tasks.

In a closely-related work, La Torre et al. [START_REF] La Torre | Context-bounded analysis of concurrent queue systems[END_REF] propose a "context-bounded" analysis of shared-memory multi-pushdown systems communicating with messagequeues. According to this approach, one "context" involves a single process reading from its queue, and posting to the queues of other processes, and the number of contexts per execution is bounded. Our work can be seen as an extension in a few ways. First, and most trivially, in their setting a process cannot post to its own message queue; this implies that at least 2k contexts must be used to simulate k phases of a single-processor program. Second, there are families of 1-phase executions which require an unbounded number of task-contexts to capture; the execution order D 1 D 2 D 3 . . . D 2n of Figure 5b is such an example. We conjecture that bounded phase and delay captures context-bounding-i.e., there exists a polynomial function f : N → N such that every k-context bounded execution of any program P is also a f (k)-phase and delay bounded execution. Finally, though phase-bounding leads to a convenient sequential encoding, we are unaware whether a similar encoding is possible for context-bounding.

Boigelot and Godefroid [START_REF] Boigelot | Symbolic verification of communication protocols with infinite state spaces using QDDs[END_REF] and Bouajjani and Habermehl [START_REF] Bouajjani | Symbolic reachability analysis of fifochannel systems with nonregular sets of configurations[END_REF] have proposed analyses of message-passing programs by computing explicit finite symbolic representations of message-queues. As our sequentialization does not represent queues explicitly, we do not restrict the content of queues to convenientlyrepresentable descriptions. Furthermore, reduction to sequential program analyses is easily implementable, and allows us to leverage highly-developed and optimized program analysis tools.

Conclusion

By introducing a novel phase-based characterization of message-passing program executions, we enable bounded program exploration which is not limited by message-queue capacity nor the number of processors. We show that the resulting phase-bounded analysis problems can be solved by concise reduction to sequential program analysis. Preliminary evidence suggests our approach is at worst competitive with known task-order respecting bounded analysis techniques, and can easily scale where those techniques quickly explode.

B Sequential Program Semantics

For expressions without program variables, we assume the existence of an evaluation function • e : Exprs → ℘(Vals) such that e = Vals. For convenience, given a processor configuration κ = g, w, q and w = , s w , we define global assignment,

κ(l ← ) def = g, , s w , q local assignment, κ • f def = g, f • w, q
append stack frame.

To further reduce clutter in the operational program semantics, we introduce a notion of context. A statement context S is a term derived from the grammar S ::= | S; s, where s ∈ Stmts. We write S[s] for the statement obtained by substituting a statement s for the unique occurrence of in S. Intuitively, a context filled with s, e.g., S[s], indicates that s is the next statement to execute in the statement sequence S[s]. Similarly, a processor configuration context C = g, , S w, q is a processor configuration whose top-most frame's statement is replaced with a statement context, and we write C[s] to denote the processor configuration g, , S[s] w, q . When e is an expression, we abbreviate e(C [skip]) by e(C).

Figure 10 defines the transition relation → S for the standard sequential program statements. The Skip rule simply steps past the skip statement. The Assume rule proceeds only when the given expression e evaluates to true. The Assign statement stores the value of a given expression in either the local variable l or the global variable g. The If-Then and If-Else rules proceeds to either the then or else branch, depending on the current valuation of the given expression e. Similarly, the Loop-Do and Loop-End rules proceed to (re-)enter the loop when the given expression e evaluates to true, and step past the loop when e evaluates to false. More interestingly, the Call rule creates a new procedure frame f by evaluating the given argument e, and places f at the top of the procedure-frame stack. The Return rule removes the top-most procedure frame from the stack, and substitutes the valuation of the return expression e into the assignment x := left below by the matching call statement. Note that the transition relation → S is non-deterministic, since the evaluation of an expression e can result in an arbitrary set of possible values.

C Full Listing of Example Programs of Section 4

The first program P 1 (i), parameterized by i ∈ N, has a single Boolean global variable b, i procedures named p 1 , . . . , p i , which assert b to be false and set b to true, and i named q 1 , . . . , q i which set b to false. Initially, P 1 (i) sets b to false, and enters a loop in which each iteration posts some p j followed by some q j . Since a q j task must be executed between each p j task, each of the assertions are guaranteed to hold. Figure 6a compares the time required to verify P 1 (i) (using the Boogie verification engine [START_REF] Barnett | Boogie: An intermediate verification language[END_REF]) for various values of i, and various bounds n on loop unrolling.Note that although every execution of P 1 (i) has only 2 phases, to explore all n loop iterations in any given execution, the size of queues must be at least 2n, since two tasks are posted per iteration.

Call v ∈ e(C) f = v, sp C[call x := p e] S -→ C[x := ] • f Return f = , S[return e] v ∈ e(C • f ) C[x := ] • f S -→ C[x := v]
Our second program P 2 is biased to support the bounded task-queue encoding. P 2 again has a single Boolean global variable b, and two procedures: p 1 asserts b to be false, sets b to true, and posts p 2 , while p 2 sets b to false and posts p 1 . Initially, the program P 2 sets b to false and posts a single p 1 task. Again here, since a p 2 task must execute between each p 1 task, each of the assertions are guaranteed to hold. Figure 6b compares the time required to verify P 2 for various bounds n on the number of tasks explored. Note that although every execution of P 2 uses only size 1 task-queues, to explore all n tasks in any given execution, the number of phases must be at least n, since each task must execute in its own phase.

D The Multi-Processor Breadth-First Scheduler

Here we define a deterministic, non-blocking, delay-accessible delaying scheduler M bfs which though perhaps odd from an operational point of view, has a very useful application: given a multi-processor message-passing program P , the phaseand delay-bounded executions of P according to M bfs are simulated by executions of a sequential program P ; furthermore, P is obtained by a simple code-to-code translation of P which does not explicitly represent pending-task queues.

Let U be a set of identifiers uniquely identifying each task along an execution with a single initially-pending task u 0 ∈ U . Our scheduler keeps a monotonically increasing phase number i ∈ N, along with an ordered task-posting tree T over nodes U , a completion-labeling √ : U → B, and a phase-labeling Φ : U → N. Initially the tree contains a single node u 0 , with Φ(u 0 ) = 0 and √ (u 0 ) = false. As additional tasks are posted, we add them as children of the posting task, in the order they are posted. Normally, the scheduler allows tasks to execute to completion; when a task does complete, the scheduler marks it as completed. When choosing the next task to execute, our scheduler selects the smallest-in breadth-first order over the task-posting tree-unexecuted task in the current phase; if there are no non-completed tasks in the current phase, the scheduler moves to the next phase. In this way, the scheduler executes all tasks in phase order, and same-phase tasks in breadth-first order of the task-posting tree.

To implement delaying, our scheduler also keeps a phase-delay counter ∆(ρ) : N for each processor ρ. Supposing an executing task u has phase-i on a processor whose phase-delay counter has current value j, the task u is treated as though it is in phase i + j. When a processor is delayed, its phase-delay counter is simply incremented; the effect is to shift all following tasks on the given processor one additional phase later. Delaying causes the currently executing task to be interrupted and resumed in the following phase.

Formally the Breadth-First Scheduler M bfs = D, empty, enabled, step, delay is defined over scheduler objects m = i, T, √ , Φ, ∆ ∈ D as described above; the initial object is empty = 0, T 0 , √ 0 , Φ 0 , ∆ 0 , where T 0 is the single-node tree with root u 0 , √ 0 (u) = false and Φ 0 (u) = 0 for all u ∈ U , and ∆ 0 (ρ) = 0 for all ρ ∈ Pids. The enabled( i, T, √ , Φ, ∆ , ξ) operation uniquely returns the processor identifier ρ of the smallest task u-according to the breadth-first order of T -such that Φ(u) + ∆(ρ) = i. The step( i 1 , T 1 , √ 1 , Φ 1 , ∆ , ξ 1 , ξ 2 ) operation for a transition τ = ξ 1 (ρ) → S ξ 2 (ρ) returns i 2 , T 2 , √ 2 , Φ 2 , ∆ such that -If τ is a Complete-step of task u, then

√ 2 = √ 1 (u → true); otherwise √ 2 = √ 1 .
-If τ is a Post-or Self-Post-step of task u posting task u , then T 2 is obtained from T 1 by adding to u new a rightmost child u , and Φ 2 = Φ 1 (u → Φ 1 (u) + ∆(u)); otherwise, T 2 = T 1 and Φ 2 = Φ 1 .

-If there no longer exists a non-completed task u on some processor ρ such that Φ 2 (u) + ∆(ρ ) = i 1 then i 2 = i 1 + 1; otherwise i 2 = i 1 .

The delay( i 1 , T, √ , Φ, ∆ 1 , ρ, ξ) operation returns i 2 , T, √ , Φ, ∆ 2 such that -∆ 2 = ∆ 1 (ρ → ∆ 1 (ρ) + 1) increments ∆ 1 's mapping for processor ρ.

-If there no longer exists a non-completed task u on some processor ρ such that Φ(u) + ∆ 2 (ρ ) = i 1 , then i 2 = i 1 + 1; otherwise i 2 = i 1 .

According to our definition, M bfs repeatedly picks a unique processor ρ to execute such that ρ is non-idle or has pending tasks, and ρ's first non-idle or pending task u has the lowest offsetted phase Φ(u) + ∆(ρ) of any task on any processor. Note that M bfs is a deterministic delaying scheduler which executes all tasks of a given phase before any task of a subsequent phase. Since M bfs must pick an enabled task so long as there are pending tasks on some processor, M bfs is non-blocking. Finally, since for any i = Φ(u) + ∆(ρ), repeatedly delaying every other processor ρ = ρ eventually increments ∆(ρ ) such that for any pending u on ρ , Φ(u ) + ∆(ρ ) > i, M bfs is delay accessible.

On non-delaying executions, M bfs essentially performs a phase-by-phase breadth-first traversal of the task-posting tree T -a tree which includes tasks across all processors. Interestingly and essentially for our sequential encoding of M bfs in Section 5, on a per-phase basis, with respect to any individual processor, the breadth-first traversal of the task-posting tree is identical to depth-first traversal. This follows from the fact that no task may have a same-processor ancestor in the same phase, and that processors do not share memory.

E Proofs to Selected Theorems

Theorem 3. The k-phase bounded state-reachability problem is undecidable for multi-processor programs with finite data domains.

Proof. We proceed by reduction from Post's correspondence problem [START_REF] Post | A variant of a recursively unsolvable problem[END_REF]: given words α 1 , . . . , α n , β 1 , . . . , β n ∈ Σ * of a finite alphabet Σ such that |Σ| ≥ 2, find

P

  Fig.1. The grammar of asynchronous message-passing programs P . Here T is an unspecified type, and e, p, and ρ range, resp., over expressions, procedure names, and processor identifiers.

Fig. 3 .Fig. 4 .

 34 Fig. 3. The multi-processor transition relation → M parameterized by a scheduler M = D, empty, enabled, step .

  Formally, a scheduler M = D, empty, enabled, step consists of a data type D of scheduler objects m ∈ D, a scheduler constructor empty ∈ D, a scheduler decision function enabled : (D × (Pids → Pconfigs)) → ℘(Pids), and a scheduler update function step : (D × (Pids → Pconfigs) × (Pids → Pconfigs)) → D. The arguments to enabled allow a scheduler to decide which processors are enabled depending on the execution history. A scheduler is deterministic when |enabled(m, ξ)| ≤ 1 for all m ∈ D and ξ : Pids → Pconfigs, and is non-blocking when for all m and ξ, if

Fig. 5 .

 5 Fig. 5. Phase-bounded executions with processors A, B, C, and D; each task's label (e.g., Ai) indicates the processor it executes on (e.g., A). Arrows indicate the posting relation, indices indicate execution order on a given processor, and dotted lines indicate phase boundaries.

Fig. 6 .

 6 Fig. 6. Time required to verify (a) the program P1(i), and (b) the program P2 with the Boogie verification engine using various encodings (bounded queues, bounded phase), and various loop unrolling bounds. Time-out is set to 100s.

Fig. 7 . 11 Fig. 8 .Fig. 9 .

 71189 Fig. 7. The delay operation.

  g/g, /l] e to evaluate the expression e in a processor configuration κ (alternatively, in a global valuation g and local valuation ) by substituting the current values for variables g and l. As these are the only program variables, the substituted expression e[g/g, /l] has no free variables. Additionally we define κ(g ← g ) def = g , w, q

  [skip] (x ← v) If-Then true ∈ e(C) C[if e then s1 else s2] S -→ C[s1] If-Else false ∈ e(C) C[if e then s1 else s2] S -→ C[s2] Loop-Do true ∈ e(C) C[while e do s] S -→ C[s; while e do s] Loop-End false ∈ e(C) C[while e do s] S -→ C[skip]

Fig. 10 .

 10 Fig. 10. The single-processor transitions relation → S for the standard sequential program statements.

  var b: bool // for j = 1, . . . , i proc pj () assert !b; b := true; return // for j = 1, . . . , i proc qj () b := false; return proc main () b := false; while do if then post p1 () else if then post p2 () .. else post pi (); if then post q1 () else if then post q2 () .. else post qi () return

In the presence of the assume statement, only the values reached in completed executions are guaranteed to be valid.

[START_REF] Abdulla | Verifying programs with unreliable channels[END_REF] We assume each task in a given execution has implicitly a unique task-identifier.

Acknowledgments

We thank Constantin Enea, Cezara Dragoi, Pierre Ganty, and the anonymous reviewers for helpful feedback.

* Partially supported by the project ANR-09-SEGI-016 Veridyc. † Supported by a Fondation Sciences Mathématiques de Paris post-doctoral fellowship.

A Syntactic Extensions Used in Our Code Translations

The following syntactic extensions are reducible to the original program syntax of Section 2.1. Here we freely assume the existence of various type-and expressionconstructors. This does not present a problem since our program semantics does not restrict the language of types nor expressions.

Multiple types. Multiple type labels T 1 , . . . , T j can be encoded by systematically replacing each T i with the sum-type T = j i=1 T i . This allows local and global variables with distinct types.

Multiple variables. Additional variables x 1 : T 1 , ..., x j : T j can be encoded with a single record-typed variable x: T , where T is the record type { f 1 : T 1 , ..., f j : T j } and all occurrences of x i are replaced by x.f i . When combined with the extension allowing multiple types, this allows each procedure to declare any number and type of local variable parameters, distinct from the number and type of global variables.

Local variable declarations. Additional (non-parameter) local variable declarations var l': T to a procedure p can be encoded by adding l' to the list of parameters, and systematically adding an initialization expression (e.g., the choice expression , or false) to the corresponding position in the list of arguments at each call site of p to ensure that l' begins correctly (un)initialized.

Unused values. Call assignments call x := p e, where x is not subsequently used, can be written as call _ := p e, where _: T is an additional unread local variable, or simpler yet as call p e.

Unused branches. if e then s else skip is abbreviated by if e then s.

Increment. Increment operations x++ are encoded as

Let bindings. Let bindings of the form let x: T = e in can be encoded by declaring x as a local variable var x: T immediately followed by an assignment x := e. This construct is used to explicate that the value of x remains constant once initialized. The binding let x: T in is encoded by the binding let x: T = in where is the choice expression.

Arrays. Finite arrays with j elements of type T can be encoded as records of type { f 1 : T , ..., f j : T }, where f 1 . . . f j are fresh names. Occurrences of terms a[i] are replaced by a.f i , and array-expressions [ e 1 , ..., e j ] are replaced by record-expressions { f 1 = e 1 , ..., f j = e j }.

For this problem instance we build the following finite-data asynchronous message-passing program with four processors ρ 0 , ρ 1 , ρ 2 , and ρ 3 : Initially, the main procedure is pending on processor ρ 0 .

In each loop iteration, main chooses a branch corresponding to an index i ∈ {1 . . . n} and posts each symbol of α i individually and in order to ρ 1 , and each symbol of β i individually and in order to ρ 2 . In this way, main sends to ρ 1 the sequence α i1 . . . α i k , and to ρ 2 the sequence β i1 . . . β i k in k loop iterations, each terminated by a last message. Each instance of the hold tasks which execute on ρ 1 and ρ 2 simply propagate their symbol to ρ 3 . Using the global variable turn, ρ 3 ensures that he only sees symbols sent from ρ 1 and ρ 2 in alternating order, starting with ρ 1 . Using the global variable prev, ρ 3 ensures that each symbol of β i1 . . . β i k sent from ρ 2 matches the previous symbol of α i1 . . . α i k sent from ρ 1 . Finally, if ρ 3 receives both terminating last' messages from ρ 1 and ρ 2 before another L-turn, then he has successfully checked the equality of sequences α i1 . . . α i k = β i1 . . . β i k . Thus, when down[L] and done[R] are both set to true, this means main was able to guess a solution i 1 . . . i k to the correspondence problem.


Hocine Menana
email: hocine.menana@univ-nantes.fr

Mouloud Féliachi
email: mouloud.feliachi@univ-nantes.fr

Electromagnetic Characterization of the CFRPs Anisotropic Conductivity: Modeling and Measurements

Keywords: 82, 20, Wt Computational modeling; simulation -81, 70, Ex Nondestructive testing: electromagnetic testing, eddy-current testing

This work deals with the characterization of the conductivity tensor of a carbon fiber reinforced polymer composite (CFRP) thin plate. We propose a contactless method based on the eddy current non destructive testing technique. The used eddy current sensor consists of a ferrite torus on which a winding is wound. The torus is of a rectangular section and contains a thin air-gap in which the thin CFRP plate is inserted. We developed analytical relations giving the longitudinal and transversal conductivities of the CFRP plate as functions of the impedances variations of the eddy current sensor, corresponding to the orientations of the carbon fibers parallel and transverse to the direction of the torus width which is much greater than its thickness. The analytical relations are developed by inverting interpolation functions of curves giving the variations of the sensor impedances as functions of the longitudinal and transversal conductivities of the CFRP plate. These curves are obtained by a numerical model based on a simplified integro-differential formulation in terms of the electric vector potential in the CFRP plate, coupled to the magnetic circuit equations in the ferrite torus. The modeling results are supported by measurements.

Introduction

Carbon fiber reinforced composites (CFRPs) present superior characteristics compared to metals; they are both robust and light weight which make them widely used especially in aerospace applications. Knowing the electrical conductivity of CFRPs, as well as the other electromagnetic properties, is of great importance for various applications such as electromagnetic compatibility [1], eddy current non-destructive testing [2], induction heating [3], and for some electromechanical applications [4]. Unfortunately, these proprieties are not known in advance; they are determined either by using the homogenization techniques [3], or by appropriate measurements once the material is constructed [5]- [6]. The four points bridge technique has been used to determine the electrical conductivity of CFRPs [5]; however, the high anisotropy and low dimensionality of these materials make this technique in many cases difficult and unreliable [6]. Contactless methods, based on resonant circuits [6] and on the physics of the eddy currents [7], have been proposed to determine the global electrical conductivity of CFRPs; however, the anisotropy has been weakly treated. In this work, an inverse model based on the physics of eddy current is developed to determine the conductivity tensor of a CFRP thin plate. The used eddy current sensor consists of a ferrite torus on which a winding is wound. The torus is of a rectangular section and contains a thin air-gap in which the thin CFRP plate is inserted. The use of the ferrite torus allows us to eliminate the perturbations due to the lift-off variations, and to impose a normal magnetic source field in the thin CFRP plate.

The experimental device and method principle

A schematic view of the experimental device is described in Figure 1. A CFRP ply of thickness δ and of unknown conductivity is positioned in the air-gap (e) of a ferrite torus of a rectangular section and on which a winding is wound. The torus winding is connected to a high precision impedance analyzer. The inner and outer radii of the torus are denoted respectively by R i and R e . The width (τ) of the torus is much greater than its thickness (R e -R i ). The numerical values of the system parameters are given in Table 1. A CFRP ply is of about 125µm thickness. It is constituted of carbon fibers of about 10µm diameter embedded in an electrically non conductive matrix. The conductivities vary between 5×10 3 and 5×10 4 S/m in the direction of the fibers, and between 10 and 10 2 S/m in the transverse direction, depending on the type and the volume fraction of the fibers [START_REF] Grimberg | Proceedings of the 8 th International Conference of the Slovenian Society for Non-Destructive Testing[END_REF]. The conductivities in the directions parallel and transverse to fibers are denoted respectively by (σ L ) and (σ T ). The insertion of the CFRP ply in the torus air-gap leads to a variation of the inductor impedance due to induced currents in the ply. As the torus is of a rectangular section and the fibers are unidirectional in the ply, this impedance variation depends on the orientation (θ) of the CFRP ply (i.e., carbon fibers) in the torus air-gap. According to the physics of the eddy currents, the longitudinal or transversal conductivity of the CFRP ply can not act alone on the impedance variation of the inductor; the cancellation of one of them implies the cancellation of the eddy currents. We propose thus to express the variations of the inductor resistance as functions of the summation of product of the longitudinal and transversal conductivities weighted by coefficients to be determined.

∑ = ∆ i i T i L i a R θ β θ α θ θ σ σ (1)
Notice that the relation (1) applies to any inductor that realizes a distinction between the directions parallel and transverse to fibers. Limited to one term, and for the orientations (θ=0° and θ=90°), the relation (1) becomes

0 0 0 0 β α σ σ T L a R ≈ ∆ , ( 1
β α β α β σ - -                   ∆         ∆ = a R a R L ( 2)
β α β α β β α α β σ - -                     ∆         ∆         ∆ = a R a R a R T (3)
The parameters a 0 , a 90 , α 0 , α 90 , β 0 and β 90 are determined from (1-a) and (1-b), by using curves giving the variations of |∆R 0 | and | ∆R 90 | as functions of σ L and σ T which are constructed by numerical simulation (figure 4). At this end, the relations (1-a) and (1-b) are linearized as follows

0 0 0 0 ln ) ln( ) ln( ) ( ln R a T L ∆ ≈ + + σ β σ α , (1-c) 90 90 90 90 ln ) ln( ) ln( ) ( ln R a T L ∆ ≈ + + σ β σ α . (1-d)
The parameters a 0 , α 0 and β 0 , as well as a 90 , α 90 and β 90 , are determined by a simple identification on three points; i.e.

The electromagnetic modeling

The evaluation of the inductor resistance variation passes through the determination of the eddy currents in the CFRP ply which is done numerically. The numerical modeling of the eddy current in the CFRP ply is coupled to the magnetic circuit equations of the inductor (torus + winding) to avoid integrating the latter in the numerical modeling. The modeled system is represented schematically in the figure 2. The influence of the CFRP ply is represented by a fictitious loop carrying a current I representing the eddy currents in the CFRP ply. The fictitious current I is calculated numerically. We assume that the magnetic flux density B is homogeneously distributed in the torus section; thus it can be expressed as follows

{ } ) ( ) 1 ( ) ( ~0 B I I N e R R B s r i e r - - + + = µ π µ µ (4)
In ( 4), µ 0 and µ r are respectively the free space permeability and the relative permeability of the ferrite torus; N is the number of turns of the torus winding; R i and R e are respectively the inner and outer radii of the torus; e is the air-gap thickness; Ĩ s is the source current and Ĩ is the current in the fictitious loop, representing the eddy current in the CFRP ply which in turn depends on the flux density. The tilde sign means that the values are complex.

The electromagnetic problem in the CFRP ply is formulated by an integro-differential equation involving the electric vector potential [START_REF] Menana | [END_REF]. The frequencies used are such that the CFRP ply thickness is much smaller than the skin depth corresponding to its maximum conductivity; in addition, the flux density is normal to its plane. The thin shell approximation is thus applied to the formulation which is reduced to involve only the normal component of the electric vector potential.

Considering the Cartesian coordinates where the CFRP ply thickness is oriented in the z axis direction, the electromagnetic problem in the CFRP ply is formulated as follows ( )

) ( ' ) ' )( ' ( ) ' )( ' ( 4 ) ( ) ( ) ( ) ( ) ( 3 ' ' 0 r B i ds r r y y r T x x r T i r T ply S z x z y z x y yy yx xy xx x y ω π δ ωµ θ ρ θ ρ θ ρ θ ρ - = - - ∂ + - ∂ +                   ∂ - ∂         ∂ - ∂ - ∫ r r (5)
In (5), S ply denote the CFRP ply surface and δ is its thickness, ω is the angular frequency, ∂ u denote the derivative with respect to the coordinate u (u=x,y), r and r' denote points identified respectively by the coordinates (x,y) and (x',y') and ' r r r ris the distance between them; T z represents the z component of the electric vector potential, i is the imaginary unit (i 2 =-1) and ρ uv (θ) (u=x,y; v=x,y) is a component of the resistivity tensor of the CFRP plate depending on the orientation (θ) of the fibers in the XY plane. The components of the eddy current density in the CFRP ply are calculated as follows

        ∂ - ∂ =         z z y x T x T y J J ~ (6)
The current in the fictitious loop, representing the effect of the eddy currents in the CFRP ply of length L, width W and thickness δ, is calculated as follows (see figure 3)

∫ ∫ = = = = 2 / 0 2 / 0 ) 2 / ( ) 2 / ( ~L y W x W y at dx J L x at dy J I δ δ (7)
The equality between the two integrals in (7), evaluated respectively at the positions x=L/2 and y=W/2, is due to the total current conservation in the plate. The equations ( 4) to (7) are to be solved iteratively. We start by evaluating (4) considering only the source current Ĩ s and then solve (5), evaluate ( 6) and ( 7) and then reevaluate (4) by introducing the updated value of the current Ĩ. The iteration process is repeated until convergence; the convergence test can be realized either on the magnetic flux density or on the electric vector potential. Once the eddy currents are determined in the CFRP ply, the variation of the inductor resistance is evaluated by the following relation, where * [ ]

∫ + + + = ∆ ply S y yy x yx y y xy x xx x s ds J J J J J J I R ) ) ( ) ( ( ) ) ( ) ( ( ) (
The system is simulated at a frequency of 100 KHz. The differential and integral terms of ( 5) are discretized by using respectively the finite difference method and the moment method. The resulting algebraic system is solved by using the bi-conjugate gradient method. Figure 3 represents zooms on the eddy currents repartitions in the CFRP ply corresponding to the orientations of the fibers in the directions parallel (θ=0) and transverse (θ=π/2) to the width τ of the torus. The influence of the anisotropy on the eddy currents repartition in the CFRP ply is noticeable; the eddy currents are stretched in the direction of the fibers (axis of the high conductivity). Figure 4 represents the variations of |∆R 0 | and | ∆R 90 | as functions of σ L and σ T in the intervals corresponding to the standard conductivities of CFRPs. We can remark that both |∆R 0 | and | ∆R 90 | vary in a significant way with the variation of the transverse conductivity σ T ; however, for a fixed σ T , when σ L >σ T , the variation of σ L leads to weak variations of |∆R 0 | and | ∆R 90 |. These results are used to identify the parameter a 0 , a 90 , α 0 , α 90 , β 0 and β 90 by using (1-c) and (1-d). The resistance variations |∆R 0 | and | ∆R 90 | by have been measured by means of the impedance analyzer at a frequency of 100 KHz, and then the conductivities σ L and σ T have been identified by using the relation ( 2) and (3). As expected, the weak variations of |∆R 0 | and | ∆R 90 | according to σ L lead to a degenerate identification of the latter by means of (2). The conductivity σ L is thus determined by (1-a) once σ T is determined by means of (3). The obtained results are given in Table 2. Figure 5 represents the measured variation of the inductor impedance (R+ iX) according to frequency, without and with the CFRP ply in the torus air-gap, with the fibers oriented in the directions parallel and transverse to the torus width. We remark that the presence of the CFRP ply and the orientation of its fibers in the torus air-gap lead to a significant variation of the resonance frequency and a significant variation of the inductor resistance at these resonance frequencies. These parameters constitute, thus, interesting features for anisotropic material characterization.

To check the validity of the identified conductivities, we compared between the calculated and measured values of the resistance variation of an air-cored coil above the characterized CFRP ply for different frequencies. A schematic view of the validation device is given in Figure 6. The coil characteristics are given in Table 1, and the comparison results are presented in Figure 7. The measured and calculated results of the air-cored coil resistance variation are globally in agreement. When looking to the low frequencies (f <700 KHz), the CFRP ply conductivities seems to be overestimated since the calculated resistance variations is greater than the measured ones. The situation is inversed at high frequencies where the measured resistance variations become greater than the calculated ones; this is due to capacitive couplings between fibers in the CFRP ply which are not taken into account in the modeling.

Conclusion

We proposed a contactless method for the characterization of the anisotropic conductivity of a thin CFRP plate. A simplified modeling is proposed and measurement results are provided. The modeling needs to be improved by taking into account the capacitive coupling in the CFRPs. We have seen that the high anisotropy of CFRPs leads to a weak dependence of the impedance variation of the inductor on the highest conductivity, which results in errors in the determination of the latter. Measurements showed that the resonance frequency of the inductor, as well as the value of its resistance at this resonance frequency, vary in a significant way with the variation of the orientation of the CFRP ply fibers in the torus air-gap. This phenomenon constitutes an interesting alternative for characterizing anisotropic materials.

Figures and tables captions

  a) and (1-b) ∆R 0 and ∆R 90 are the resistance variations of the inductor corresponding respectively to the positions where the fibers are oriented in the direction of the torus width and in the transverse direction. The relations (1-a) and (1-b) are inverted to express the conductivities σ L and σ T of the CFRP plate as functions of ∆R 0 and ∆R 90 ; after some manipulations, we obtain the following relations

  solving 1-c and 1-d for three chosen values of |∆R 0 | and | ∆R 90 | taken from figure 4.

Fig. 1 .Fig. 2 .Fig. 3 .Fig. 4 .Fig. 5 .Fig. 6 .Fig. 7 .

 1234567 Fig. 1. Schematic view of the experimental device

Figures

  Figures

Fig. 3

 3 Fig. 3

Fig. 4

 4 Fig. 4

  Fig. 5

Fig. 6

 6 Fig. 6

Fig. 7

 7 Fig. 7

Table 1 .

 1 The numerical values of the system parameters

Table 2 .

 2 Measurement and identification results

Acknowledgments

Work supported by Nantes Métropole.

Tables