F Gratecap

M Girard

S Marya

G Racineux

Exploring material flow in friction stir welding: Tool eccentricity and formation of banded structures

Keywords: Friction stir welding, Eccentricity, Mixing, Material flow, Band spacing

Surface striations with spacing equal to feed rate per rotation and banded structures in the weld nugget are some of the striking features of friction stir welding. However, their formation is still subject to some debate. This study contributes to comprehend their formation by evaluating the possible role played by the eccentricity of the tool during the welding of an aluminium alloy and using a plasticine as its analogue. The eccentric movement is visualized to generate both surface and bulk striations in plasticine. By voluntarily using non optimized welding conditions on aluminium, the material flow has been deduced and confronted with direct visual observations through high speed camera on plasticine. In the non closed section of the weld, two lobes each with thickness equal to feed rate per rotation were observed. First lobe corresponds to flow induced by the pin and tests on plasticine showed constant volume displacement per rotation for a given tool eccentricity. The second lobe is generated by material flow from under the shoulder back to the rear of the pin. The assembling or not of these two lobs behind the pin can explain some of the characteristic patterns observed in the weld nugget such as onion rings, oxide dispersions and cavities.

Introduction

Friction stir welding ever since its invention has undergone extensive academic investigations and demonstrated its viability in industrial applications. Investigations on friction stir welding have so far mainly focused on developing tools and procedures for making welds in a variety of alloys, characterizing the properties of welds for allowable designs and investigating the possible flow mechanisms in and around the weld nugget. The tool design is of prime importance in explicating temperature fields and material movement patterns in FSW. The process involves a milling type machine with different degrees of freedom. The tool generally consists of a shank, a shoulder and a pin that is adjusted to a certain height inside the work pieces to be welded. The tool is rotated and a relative displacement motion along the weld line is generated. Different tool geometries from a simple pin to complex threads have been proposed and investigated. The basic idea behind the process is the generation of heat by friction between the tool and workpieces and its plasticization with rubbing from the moving tool. The induced material deformation originates by dry rubbing, lubricated and sticking depending on the processing conditions and material chemistry. Under optimized conditions, the hole left behind by the advancing tool is filled by material input generated by the forward motion of the tool. As the tool rotates and advances in a very continuous manner, the process is supposed to leave a weld surface without any print marks other than those arising from the tool's surface itself. However, one of the striking features of FSW nuggets is the presence of striations with spacing equal to forward motion (V f , mm/minute) per rotation (N, rotations/minute) [START_REF] Chen | Shear flow and formation of nugget zone during friction stir welding of aluminium alloy 5083o[END_REF]2]. Other types of striations and patterns in the weld nugget have been observed and characterized. For example, the formation of typical onion F. Gratecap : M. Girard : S. Marya : G. Racineux (*) Institut de recherche en génie civil et mécanique (GeM) -UMR CNRS 6183, Ecole Centrale de Nantes -1 rue de la Noë, 44321, Nantes cedex 3, France e-mail: guillaume.racineux@ec-nantes.fr rings in the transverse cross-section of the nugget, the formation of striations in the longitudinal sections and to a lesser extent, the distribution of typical oxide patterns in different sections of the weld [START_REF] Zettler | The effect of process parameters and tool geometry in friction stir welding of the alloys az31 and az61[END_REF][START_REF] Chen | The investigation of typical welding defects for 5456 aluminum alloy friction stir welds[END_REF]. The contribution of the material flow patterns around the moving tool on the formation of different banded structures of the weld nugget has been proposed. Investigations of flow patterns involve deductive conclusions from post mortem specimens by materializing tracer distribution in the weld or by tool piece frozen experiments [START_REF] Dickerson | A weld marker technique for flow visualization in friction stir welding[END_REF]. Though tracer techniques suffer from possible discontinuities generated in material flow by their presence, they still offer the best possible channel to study material flow. For example, Colligan et al. [START_REF] Colligan | Material flow behaviour durinf friction stir welding of aluminium[END_REF] used the steel shot tracer technique in conjunction with a screw threaded pin and concluded limited stirring of the surface material and extrusion of the other materials on the retreating side of the nugget. Lorrain et al. [START_REF] Lorrain | Understanding the material flow path of friction stir welding process using unthreaded tools[END_REF] made welds using two different tool profiles and copper foils as material marker. Zettler et al. [2] used Ti powder markers for different aluminium alloys and investigated material flow by using different tool geometries. They showed a good correlation between the feed rate per rotation and band spacing whatever the tool geometry, though significant temperature field differences could be found between different tool geometries. In fact, the published literature suggests that striations arise from separation of the material flow pulled in to the stir zone from either side of the joint. From tracer techniques, the flow patterns are visualized as shearing from the advancing side followed by the material flow from the retreating side in far field of the pin. Thus, in the nearest flow field, closest to the pin, the material originates from the advancing side and is forced to flow to the other side of the joint. Material from the retreating side enters the stir zone as a separate entity from behind the tool pin due to the lower pressure generated there by the advancing tool. The two flows do not mix in real sense of term and generate layered structures. In case of dissimilar joints Larsson et al. [START_REF] Larsson | Joining of dissimilar al-alloys by friction stir welding[END_REF], while intimating bonding is observed, the analysis of composition profiles generally shows only a very limited mixing, which supports the flow observations from tracer studies. Further onion rings, striking features of FSW seem to result from extrusion of the cylindrical sheets of the material during each rotation of the tool and according to Krishnan [START_REF] Krishnan | On the formation of onion rings in friction stir welding[END_REF] explains one to one relation between the spacing of the marking and forward motion of the tool per rotation. The basic ideas derived from observed or deduced flow patterns support stirring and extrusion as dominating mechanisms. Regarding the oxide patterns reported by Jene et al. [START_REF] Jene | Oxide fragments in friction stir welds-distribution and effects on crack initiation[END_REF], a good correlation between the oscillating oxides and the measured forces in transverse and longitudinal sections of wrought Al alloy AA5454 has been reported. The periodic length of the measured forces and the oxide patterns correlated closely with the feed rate per rotation (FPR). However the cause and effect relation is not clearly established. It is not sure if in the absence of the oxide patterns, the oscillating forces would vanish or not. The force plot suggested a circular motion of the tool, caused by the unbalanced mass of the tool or the asymmetric friction stir itself due to difference on the advancing and retreating sides. The effect of eccentric weld tool, whatever the reason, on the formation of different pattern characteristics remains comparatively scarce. This work proposes to contribute to this problematic by reporting a number of observations made on FSW welds on aluminium substantiated by investigations on plasticine as analogue material. It is worth mentioning that plasticine has previously been used to model forging and extrusion characteristics and offers a simple way to understand possible flow characteristics of metals. More recently Leichty et al. [11] used plasticine to explore flow characteristics in friction stir welding as it shows typical stress strain behaviour with temperature and strain rates as metals. Further due to very low stresses encountered by the tool in plasticine, any observed eccentricity would arise mainly from tooling itself. In case of metals, reaction forces between metal and tool may account for some generated eccentricity. Recently, Kumar et al. [START_REF] Kumar | The role of friction stir welding tool on materialflow and weld formation[END_REF] proposed a special experiment, where the interaction of the friction stir welding tool with the base material is continuously increased. The results show that there are two different modes of material flow regimes involved in the friction stir weld formation; namely "pindriven flow" and "shoulder-driven flow". The objective of this study is thus to report the effect of tool eccentricity on observed flow patterns in an aluminium alloy and get better insight of their formation by matching with those on plasticine in carefully designed set of experiments. The reported results point out a very strong similarity of observed patterns and direct observation of flow in plasticine gives interesting qualitative clues to possible material flow in aluminium during the friction stir welding process.

Experimental method

The work is mainly focused on 4 mm thick aluminium alloy 2017-T4 plates and plasticine of similar thickness obtained by rolling. In this study four tools were used. The first one is a conventional FSW tool made from quenched & tempered steel with a 10 mm radius shoulder and 3 mm radius conical smooth pin (Fig. 1a). The second one is a cylindrical tool without any shoulder (Fig. 1b). The third and fourth was cylindrical with 3 and 4 milled faces at the bottom (Fig. 1c andd). Preliminary studies showed strong similarities in the profiles (surface markings, end hole.) of stirred zones between FSW of aluminium and plasticine for a given tool geometry but with understandably different working conditions. The FSW experiments were carried out on a GSP numerical milling machine. The FSW experiments are displacement controlled. Aluminium plates are rigidly fixed between clamped bars and care is taken to ensure that no lateral movement during work occurs. In order to minimize heat loss through the machine table, the workpieces were clamped on a backing plate in titanium. The heat loss effect is insignificant for plasticine, but becomes measurable for aluminium with high thermal conductivity. Gratecap [START_REF] Kumar | The role of friction stir welding tool on materialflow and weld formation[END_REF] reported that the presence of titanium supporting plate improved the consistency of test results on aluminium plates and this set up was maintained as a conservative measure in subsequent investigations. However, as the objective of this study is to understand material flow during FSW process by using plasticine, the presence of backing Ti plate would not affect the results. Temperature rise if any is quite negligible during tool advance in plasticine in contrast to aluminium. The FSW tool is fixed in the vertical broche of the milling machine with screws threaded across the thickness of the broche. The striations and flow patterns are observed by using optical microscopy and scanning electron microscope for aluminium. A direct observation of flow patterns in plasticine is looked by high speed camera with capture rate of 125 pictures/s. In general, the advance rate V f ranged from 70 to 200 mm/min and rotational speeds N was varied from 400 to 850 rpm.

Results and discussion

Typical striations of Fig. 2a on aluminium with spacing equal to feed per rotation (V f /N=0.44 mm/round) are observed and this is the general case with all experiments with feed rates and rotational speeds varying respectively from 50 to 300 mm/min and from 200 to 1000 rpm. Even if the weld parameters are slightly varied from the optimized values, these striations always prevail, but the weld may show some continuous unfilled open tunnel in the terminal part of the weld. The bottom as well as the tunnel side walls (Fig. 2b andc) shows striations with spacing of the feed per For the introduction of the mentioned tunnel defect in the stir zone, the tool was slightly lifted above the workpiece from its optimized value that was determined from tool-processing condition window [START_REF] Gratecap | Contribution to friction stir welding process[END_REF][START_REF] Gratecap | A simple methodology to define conical tool geometry and welding parameters in friction stir welding[END_REF]. By slightly lifting the tool above the workpiece, the forging conditions become insufficient to completely fill the rear of the weld and an open tunnel close to the advancing side was formed.

Further, the left terminal hole, once the tool is removed, depicts two weld lobes (Fig. 3a) with typical striations. The lobes seem to form by the addition or superposition of layers of thickness representing the feed per rotation (Fig. 3c andd). The observation of the lobes is first of its kind and certainly gives clue to what happens during the friction stir welding. As shown later on, the lower and the upper lobe are respectively generated by flow around the pin and the upper shoulder which transfers the material from front to the rear of the advancing tool. Figure 3b shows a very similar profile in the hole left by tool removal in the plasticine. Here, the conditions are different from those used for aluminium. For a better inside on these striations and lobes, simulated weld runs were made on a bloc of plasticine. In one series of the experiments, only pin was used, while in other series the whole FSW tool containing both pin and shoulder were used so as to differentiate the role of pin and shoulder. Figure 4 shows four stills taken from high speed camera. In plasticine, without shoulder, the rotation and translation of the tool yields an opening of width of size of the pin diameter and on the retreating side, a succession of layers similar to striations can be seen. The observations of the tool rotation and stacking of layers by using high speed camera revealed that the pin did not execute a completely concentric rotation, but showed an eccentricity of 0,03 mm. This is small, but its effect could be perceived as shown in Fig. 4 where after the formation of a layer, follows a period in which tool looses contact with the formed layer i.e., a gap is generated and then again a new layer is built up on the previous one. A succession of layer additions in a discontinuous manner generates the typical striation pattern on the side walls with spacing of the feed rate per rotation.

The temporary creation of gap is related to the eccentricity, a period during which the pin suddenly moves to the forward side of tool translation. During each rotation, the eccentricity of the rotating tool pin generates additional oscillatory movement in the stir zone as illustrated by the evolving position of the blue circle in Fig. 5 with reference to the red circle. The process is thus discontinuous in spite of the continuous constant advance of the tool carrying system. Here again, the number of stacked layers is governed by feed rate per rotation as shown schematically by Fig. 5. Each stacked layer can be further divided if the circular pin sections are replaced by milled facets. For example, Fig. 6 depicts three layers for three milled faces and four layers with four milled faces pin. As the number of stacked layers for a given feed rate per rotation increases with the number of milled facets, the separation between layers decreases and for a circular pin, it is expected that the number of layers would tend towards infinite which would make individual layers difficult to differentiate. This is what happens with non faceted pin where only visible layers are those generated by feed rate per rotation. Three or four striations observed in Fig. 6 occupy the same space as the layer of the circular pin without milled faces. When pin with shoulder is used to simulate the weld, then the opening left by pin motion is closed and observation of the flow patterns is impossible. The formation of lobes both in aluminium and in plasticine reported earlier (Fig. 3) supports the idea that the latter could help further in understanding material flow generated by the rotating tool. For example, it is difficult to study in-situ material flow in aluminium other than by tracer technique and post mortem metallographic analysis. On the contrary, plasticine offers alternative route to visualize material flow in FSW, though the working conditions would be completely different. For example, temperature gradients prevail in metals during FSW, but these are quasi inexistent in plasticine. The macroscopic process mechanisms were thus investigated on plasticine by designing the experimental procedure such that camera could see the motion of the tool and the material flow.

In order to visualize macroscopic phenomenon occurring in FSW, the tool pin was inserted about three quarters into the plasticine so as to leave an open space for observation. Figure 7a shows an extracted image from the film and one can notice the material coming down from the front to the back via the shoulder. Further, material layered by the pin comes from the retreating side and junction between these two lobes generates the weld. These lobes and stacking layer profiles are also seen on the side wall of the tunnel on aluminium (Fig. 7b andc). There is a very close similarity between the lobe formation seen through video on plasticine and frozen observation of the tunnel wall in aluminium. The lower lobe created by the pin is difficult to discern from the images of the Fig. 7a, but its presence is much more precisely seen from the layered structures on the retreating side when the shoulder of the tool is removed i.e., tool is a simple pin. With this arrangement, the tool leaves a open tunnel with visible side walls that contain typical stacked structures observed here.

Figure 8 shows layered structures on the retreating side as related with advance per rotation. With its reduction, the layers become more and more difficult to distinguish as can be seen when Vf/N passes from 2 to 0, 125 mm/min. Besides the layers, the Fig. 8 shows plasticine extruded well above the base on the retreating side. The presence of the shoulder blocks this material extrusion generated by the pin itself which then fills in the trench observed without shoulder. The stacks of layers are formed at each rotation of the tool which in our study revealed discontinuous feed following its eccentric rotation. The discontinuous nature of material flow results in three dimensional discontinuities which lead to different stack- ing structures revealed here on aluminium and on plasticine. From metallographic studies, Krishnan [START_REF] Krishnan | On the formation of onion rings in friction stir welding[END_REF] have attributed the onion rings observed in the cross-section of the weld nugget to such stacking layers. The formation of the surface striations is markedly explicated by the schematic diagram of Fig. 9. Due to the tilt and eccentric motion of the tool, the surface effect is as well transferred to the material in its volume which explains different banded or layered structures. Whatever the tool's eccentricity, the spacing will remain constant for a given feed per rotation. However, the quantity of material involved for a given feed per rotation will increase with increasing eccentricity. As spacing between striations depends mainly on advance per rotation, more material involvement with increasing eccentricity would generate thicker layers. The Fig. 9a also supports the up and down movement related to the tool eccentricity as striations on the bottom of the tunnel can (Fig. 7) only be generated by variable contact in between plasticine and the pin.

The stacking of layers and the formation of lobes with stacked layer structures can contribute to explicit some published results. For example, Fig. 10

Conclusions

The present study on an aluminium alloy supplemented with tests on plasticine show some interesting features that substantiate results of other researchers. One of the important features detected by high speed camera has been that the tool as a whole has an eccentric rotation movement. This generates oscillatory movement and leads to a discontinuous welding process from material perspectives in an otherwise continuous rotation translation picture of the tool. The origin of the eccentricity in the present case stem from the FSW tool attachment in the tool holder of the milling machine. In other cases, the eccentricity may also be caused by the different transverse forces from the advancing and retreating sides of the weld. This machine tool fixation aspect has very been often been neglected in explicating different flow theories around the pin and shoulder.

V f /N = 2 mm/tr V f /N = 0,5 mm/tr V f /N = 0,125 mm/tr The eccentricity reported here is very small, less than 0, 03 mm. The principal results may be summarized as follows:

1. Surface striations and stacked layer structures on different sections of the weld seem to correlate with the oscillatory motion of the rotating and advancing tool. 2. The above mentioned oscillatory motion is caused by its eccentric positioning with reference to the spindle as shown by high speed video captures of tool's advance in plasticine. 3. The eccentricity of the tool may result from the architecture of tool fixation in the broche of FSW machine alone or combined with the effect of the restraining forces when the workpiece material is resistant.

4. Spacing between stacked layers and striations is equal to advance per rotation 5. Stacking of the layers is discontinuous and correlated with eccentric forward-backward motion of the tool 6. High speed camera observations on plasticine showed two main tracks of material flow. One under the tool shoulder and the second one generated by the pin 7. These two types of material flow generates two lobes, each one with stacked structures 8. Welding is successfully accomplished when the processing conditions are such that the two lobes completely merge behind the advancing tool. 9. The two lobes, each with stacked layers, can explain the formation of reported wavy oxide patterns on aluminium. Though further studies are still required, this work underlines the importance of the machine tool alignments and fixtures in the welding process. This aspect has often been neglected in the past.

Fig. 1 aFig. 2

 12 Fig. 1 a Conical tool with β=21.7°, γ=8°, h=3.9mm, D s =20mm, D p =6mm b Cylindrical tool c tool with three milled faces d tool with four milled faces pin

Fig. 3 a

 3 Fig. 3 a Observation of two lobes (upper and lower marked by circled arrows) in the terminal hole left by tool (b) similar lobes on plasticine (c) striations on the lower (d) and upper lobe

Fig. 4

 4 Fig. 4 High speed video sequence (125 pictures/s) illustrating the formation of successive stacked layers on the simulated retreating side on plasticine. Tool 2 without shoulder (Fig. 1b)

Fig. 6 Fig. 5

 65 Fig. 6 Effect of milled faces on conical pin with the same Vf/N=1 mm/round. Each stacked layer is segmented into 3 (a) or four 4 (b) occupying the same space equal to Vf/N

Fig. 8 Fig. 7 a

 87 Fig. 8 Striations on the retreating side of plasticine with tool 2. Spacing equals Vf/N. Pin extruded plasticine can be seen beyond the thickness of the experimental plate

Fig. 9 Fig. 10 a

 910 Fig. 9 Schematics of surface striations following up and down motion of the tilted tool on plasticine. Striation thickness increases with tool eccentricity and tilt at constant Vf/N


Keywords: algèbre qualitative, bioinformatique, modèles statiques

Diagrammes de Décision pour des Modèles Biologiques Qualitatifs.

Résumé : La modélisation de grand reseaux biologiques (métaboliques, géniques ou mixtes) demande d'interpréter de nombreuses connaissances hétérogènes accumulées par les biologistes. En ce qui concerne le comportement de ces réseaux, ces informations sont le plus souvent de nature qualitative. De plus, beaucoup d'observations portent sur une comparaison entre une ou plusieurs situations expérimentales où les systèmes biologiques peuvent être considérés en équilibre. Dans plusieurs papiers précédents, nous avons montré que les variations observées entre deux situations expérimentales, pouvaient être considérées comme résultant de déplacements d'équilibre. Les signes de ces variations sont soumises à des contraintes dans l'algèbre des signes qui permettent de confronter, dans une certaine mesure, modèles qualitatifs et expérimentations. Dans ce rapport nous présentons notre implémentation des équations qualitatives comme des polynomes sur des corps finis. Des algorithmes spécifiques sont décrits en détail.

Qualitative variational models in integrative biology

Several families of mathematical models of biological network have been proposed. One of the oldest one is the family of models based on differential equations. These models are rooted in the studies of the kinetic of individual reactions in biochemistry. The differential model of a network is obtained by composing differential equations modeling individual reaction, assuming additivity of metabolite fluxes. Extensions have been proposed to gene activation and signalisation.

Models derived from differential equations

Despite the number of studies devoted to enzyme activity, only a a few reaction kinetic constant values are known. Moreover, the measurements of kinetic constants are in general made in vitro and it is not established that the values are the same in vivo. In [START_REF] Siegel | Qualitative analysis of the relation between dna microarray data and behavioral models of regulation networks[END_REF], [START_REF] Veber | Complex qualitative models in biology: A new approach[END_REF] a qualitative model was proposed which overcome to some extend, the lack of quantitative data. Let us recall briefly the derivation of such a model.

We consider a network of interacting cellular constituents. These products may be proteins, RNA transcripts or metabolites for instance. The state vector X denotes the concentration of products in the system. X is assumed to evolve according to the following differential equation:

dX dt = F(X) ( 1
)
where F is an (unknown) non linear function. A steady state X eq of the system is a solution of the algebraic equation F(X eq ) = 0

(2)

An experiment on the system consist in a comparison between two situations of interest and is modelized as an equilibrium shift. As a result to a change in the parameters (experimental settings), the system goes from its initial stable steady state to another one. In the end, we observe ,the sign of the total variation in concentration for a subset of products. Our aim here is to understand how variation signs are constrained by the interactions between products.

Let us denote s( j, i) = sgn( ∂ F i ∂ X j ) where F i is the ith component of F. We assume that s(i, i) = -, i.e. the diagonal of the Jacobian matrix of F are strictly negative. For a detailed discussion of this hypothesis see [START_REF] Siegel | Qualitative analysis of the relation between dna microarray data and behavioral models of regulation networks[END_REF]. Moreover the sign of each entry ∂ F i ∂ X j of the jacobian matrix is assumed to be constant if this entry is not null. Taking the total derivative at an equilibrium point, we get:

dX i = - ∂ F i ∂ X i -1 ∑ k∈pred(i) ∂ F i ∂ X k dX k ( 3
)
Under mild conditions discussed in [START_REF] Siegel | Qualitative analysis of the relation between dna microarray data and behavioral models of regulation networks[END_REF], it is possible to derive an equation on signs of concentration variations between two equilibria points:

s(∆X i ) ≈ ∑ k∈pred(i) s(k, i)s(∆X k ) (4)
where s(∆X k ) denotes the sign of the total variation of the concentration of kth chemical specie. Notice that these sign equations are valid not only for small variations but also for fairly large ones. These equations in the sign algebra {+, -, ?} where ? stands for an indeterminate sign, constitutes a model for qualitative variations. It is comparable with some models used in economics for comparative static analysis.

Models given in terms of influences

It often happens that biological knownledge is directly expressed in term of influences. Assertions like " when specie A concentration increases then specie B concentration decreases" are very common and some data bases such as RegulonDB( [START_REF] Salgado | RegulonDB (version 5.0): Escherichia coli K-12 transcriptional regulatory network, operon organization, and growth conditions[END_REF]) collect such type of knownledge. This type of observation is very common in gene interaction description.

It may be hypothesized that this graph results from the observation of the variations between two equilibria of a differential model. This kind of knownledge is often visualised as a graph named an interaction graph. Nodes represents molecules and edges represents relations between molecules. In this case, edges are labeled by signs and we get an signed influence interaction graph.

A major difference with qualitative models derived from differential equation, is the composition law of influences on one component (specie concentration, gene activity ...) of the model. Although a simple additive rule may be encountered, in many cases, the composition of influences targetting a component is a complex boolean function and varies from one chemical specie to another one. Think for example of a trascription factor which might induce the expression of a gene g 1 and repress the expression of another gene g 2 .

Nevertheless, gathering this kind of informations on a set of biological component ended up in a set of constraints on a set of qualitative quantities. Very often these quantities belongs to a sign algebra {+, -, ?} but in some cases several levels of concentrations are considered. The technique presented in this report applies also to these more complex qualitative algebras. For a brief presentation of sign algebra see [START_REF] Travé-Massuyès | Modèles et raisonnements qualitatifs[END_REF], [START_REF] Dormoy | Controlling qualitative resolution[END_REF].

A qualitative dynamical model

A qualitative dynamical model has been proposed for gene networks [START_REF] De | Qualitative simulation of genetic regulatory networks using piecewise-linear models[END_REF]. This model is based on the assomption that gene expressions can be modelled by a piecewise constant linear differential equation. The bounded domain of concentration vectors is partitioned with a finite number of hypercubes. On each hypercube, the coefficients of the linear model are constant. Moreover, in this equation, all the components are independant.

For each gene, the hypercube partition induces a concentration domain partition by intervals. Each interval corresponds to some level in the concentration of the gene product. So qualitative values can be associated to these levels and the evolution of the levels of gene expression is given by a family of relations, one for each gene g:

F g (X, x g ) (5)
Where X is the present vector of gene expression levels, x g is the next level of expression of gene g. It must be mentioned that given an expression level of the gene, there are several admissible next expression level for each gene. So the dynamical system is intrinsiquely non deterministic. Although this report is devoted to static models, we would like to mention that the encoding of qualitative systems on finite fields can be used for studying dynamical systems. In particular model checking can be implemented in this coding and was already used for circuit verification and real time program verification [START_REF] Marchand | Formal verification of programs specified with signal : Application to a power transformer station controller[END_REF].

Finite fields and the coding of qualitative equations

Finite fields are of widespread use in computer science. There are used in circuit verification, coding, cryptography and many other domains. Their use in circuit verification was mad possible by the discovery of an efficient representation of polynomial functions on the field / 2 [START_REF] Bryan | Graph-based algorithm for boolean function manipulation[END_REF]. This representation using a data structure known as BDD (binary decision diagrams), is based on the well known Shannon decomposition of boolean functions. In [] a similar representation of polynomial functions on / 3 = {0, 1, -1} was proposed and used for the verification of signal processing programs.

Finite fields

Finite fields or Galois fields are fields with a finite number of elements. The number of elements is necessary of the form p n where p is a prime. For a quick look at finite field theory see [START_REF]Wikipedia: Finite fields[END_REF]. A Galois field of order p n is traditionaly denoted GF(p n ).

For representing qualitative values, it is necessary to chose a Galois field with order in accordance with the finite number of values taken by each variables. In this report we will be mostly interested by simple sign algebra where each variable take value in the set {+, -, ?}. This is in relation with our biological application and all the algorithms presented here are implemented in this case. In this setting, the Galois field / 3 = GF(3 1 ) is our favorite.

Polynomial functions on finite fields

Assuming a finite field K being choosen to represent the qualitative values of interest, then every constraint, every evolution equation on the qualitative values are described by finite families of functions from K n into K. Due to the finiteness of the field, we will see soon that every such function tuns out to be a polynomial function.

The ring of polynomial functions in several variables on a finite field K of order k is itself finite. It doesn't coincide with the ring of polynomials which is infinite. It is also a finite vector space on the field K.

Let X = (x 1 , ..., x n ) a set of n variables and α an element of K. The α-Lagrange polynomial in x i is:

L α (x i ) = ∏ β ∈K;β =α (x i -β ) ∏ β ∈K;β =α (α -β ) (6)
The Lagrange polynomials satisfy:

L α (β ) = 0 if α = β and L α (α) = 1.
The family of polynomials

L (α 1 ,...,α n ) = L α 1 (x 1 ) × . . . × L α n (x n ) (7)
for all (α 1 , . . . , α n ) in K n is a basis of the vector space of polynomial functions with n variables on K.

As a consequence, it is easy to show that every function with several arguments on a finite field K is a polynomial function. If f is defined on K n with value in K, we get:

f (x 1 , . . . , x n ) = ∑ (α 1 ,...,α n )∈K n f (α 1 , . . . , α n )L α 1 (x 1 ) × . . . × L α n (x n ) ( 8
)
This equality is a direct consequence of the properties of Lagrange polynomials.

It is also possible to generalize Shannon decomposition of boolean functions to polynomial functions on finite field. Given f a (polynomial) function on K n , let us denote f | α (x 2 , . . . , x n ) the polynomial function on K n-1 resulting from the substitution of the first variable x 1 by the value α. Notice that this definition assumes an order on the variables.

It is then straightforward to check the generalisation of the Shannon decomposition:

f (x 1 , . . . , x n ) = ∑ α∈K L α (x 1 ) f | α (x 2 , . . . , x n ) (9)
If at least two polynomial functions f | α are different, x 1 is called the rootvar of f and denoted rootvar( f ). If the three functions f | α are all equal, then f is independant from x 1 . Shannon decomposition is the starting point to the representation of polynomial functions on finite fields as Direct Acyclic Graphs (DAG). The first step is to derive a tree representation from the generelized Shannon decomposition by applying it recursively. The second step follows from the remark that it often happends that many subtrees are identical. Then eliminate this redondancy by representing each subtree only once. Figure 1 illustrates the process for K = / 3 . For the original idea due to Bryant, have a look at Wikipedia [START_REF]Wikipedia: Binary decision diagrams[END_REF].

Coding of qualitative equations

Our goal is to implement efficient computations on qualitative algebra. For illustrating the method we will use the sign algebra {+, -, ?}. Application of the method to other finite qualitative algebra is similar. A sign algebra comes not only with symbols such that +, -, ? but also with operations on symbols defining a specific calculus. The sign algebra is endowed with two composition laws: qualitative sum and qualitative product.There are specified by the two tables:

⊕ + -? + + ? ? -? -? ? ? ? ?

⊗ + -? + + -? --+ ? ? ? ? ?

In addition, sign algebra has a relation, improperly named equality, denoted as ≈ and defined as follows:

≈ + -? + T F T -F T T ? T T T
It has special properties. In particular the qualitative equality is not an equivalence relation, since it is not transitive. This implies that computations in qualitative algebra must be carried with care. The two major properties to be emphasised from the above tables are:

• if a term of a sum is indeterminate (?) then the whole sum is indeterminate.

• if one hand of a qualitative equality is indeterminate, then the equality is satisfied whatever the value of the other hand is.

A solution of a qualitative system is a vector without indeterminate component satisfying the equations.

The coding We code + as 1, -as -1 and ? as 0. In order to code qualitative equations, we need to code the qualitative sum and the qualitative product. A straightforward check using value enumeration shows that the qualitative sum X ⊕Y can be coded as the polynomial function sq(X,Y ) def = XY (X +Y ) while the qualitative product X ⊗Y is simply code by the product of polynomial functions. Qualitative equations X ≈ Y are translated into polynomial equations following the pattern eq( f , g) def = f g( f -g). Solutions are obtained by decoding solutions of the polynomial equation eq( f , g) = 0.

In fact the implementation of the qualitative equations doesn't use the sum and product of polynomial functions. We take advantage of very usefull property of Shannon decomposition. Given ∆ any binary operation defined on polynomial functions, then we have:

f 1 ∆ f 2 = ∑ α∈K L α (x 1 ) f 1 | α ∆ f 2 | α (10)
So to implement a binary operation it is sufficient to have a generic implementation of binary operations and a table for the constants. Qualitative operations are implemented with the tables 2.3. With this coding, every qualitative system has a solution if and only if the corresponding polynomial system has a solution without null component. Null solutions are excluded since ? solutions are excluded for qualitative systems. In general we will have to add polynomial equations X 2 i = 1 for each variable X i , to insure this condition.

Former approach of qualitative equations

Sign algebra have allready been used in qualitative modeling applied to various plants enginered by humans. Several heuristics were proposed for the resolution of qualitative systems. For linear systems, set of rules have been designed [START_REF] Dormoy | Controlling qualitative resolution[END_REF]. This set is complete: it allows for finding every solution. It is also sound: every solution found by applying these rules is correct. The rules are based on an adaptation of Gaussian elimination. But a systematic use of elimination is not possible as it is the case for linear systems over a field. Instead, only heuristics exists for choosing the equation and the rule to apply on it. In case of dead-end when no more rule can apply, it is necessary to backtrack to the last decision made. So, programs implementing qualitative resolution were not very efficient in general and only problems of small size can be resolved in reasonable time. Our approach, inspired by circuit verification techniques, enlarge the range of application of qualitative algebra.

Problems from biological models

We have briefly presented in section 1 how static qualitative models may be used for modeling biological systems. In the present state of biological experimental technologies, it is possible to have qualitative measurements on a great number of variables. For exemples, microarray techniques allow for the simultaneous observation of the activity of thousand of genes. The challenge for bioinformaticians is to integrate numerous data from various experiments together with previous knownledge on the biological system of interest. One way to integrate biological previous knonledge is to build a static qualitative model.

Moreover, new experimentation techniques known under the generic name of hight throuput techniques became recently available in molecular biologie. They allow for simultaneous observations of many variables. With the coding of qualitative equations as polynomial functions, we have a new technique to deal with models which scale with the new technologies.

Simple problems

Consistency Given a model (a set of qualitative equation E ), the first question coming to mind is: is my model coherent?. A necessary condition for this is to have non trivial solutions. That is, in the qualitative domain, at least one solution as defined previously. A solution is obtained by defining a function encoding the conjunction of equations. In / 3 such a function is, for two equations P, Q:

f (P, Q) = (P 2 + Q 2 ) 2 (11)
since f (P, Q) = 0 if and only if P = 0 and Q = 0. This is easily proven by simple enumeration or with the help of Fermat little theorem. For any finite field a similar function is defined by a table such that f (α,

β ) = 1 if (α, β ) = (0, 0)
and f (0, 0) = 0. Applying such a function repetively replaces a set or system of equations by a unique equation

f (p 1 , f (p 2 , . . . , f (p n 1 , p n ) . . .)) (12)
with the same set of solutions. 1So checking the consistency of a static qualitative model turns out in checking wether this unique equation is the constant polynomial function 1 or not. If it is 1, there is no solution and the model is not consistent. Recall that we add the equations X 2 i = 1 to the model equations E ) to eliminate undeterminate solutions. Many practical questions are solved in the same way such that the compatibility of a model with experimental data or the efficiency of an experiment. All these questions are about the set of solutions of the model.

Prediction

One of the main role of models in science is to allow for the computation of predictions. When qualitative models are in use, the notion of prediction from the model is not quite clear. We cannot consider every solution of the system of equations as a prediction. Consider a biological system modeled as a static model on sign algebra with each variable x representing the variation of a chemical specie concentration. If the model has solutions with x = + and solutions with x = -this translates to x may increase or decrease. This can hardly be considered as a prediction!. This motivates the definition of hard components.

A component x of a vector X is a hard component for the system of qualitative equations E if for each solution (α 1 , . . . , x = α, . . . , α n ) of E , x takes the same value.

Hard components is a well known concept in qualitative modeling. It has no counterpart in polynomial function theory. So we had to develop an original algorithm to compute hard components of a polynomial function.

Experimental error

It may happen, and for a qualitative model it is the most interesting case, that experimental data do not fit the model. More precisely, the model obtained by substituing experimental data for the corresponding variables in the model E is not consistent (i.e. has no solution). In such a situation it may be either that the model is wrong, either some experimental data are questionable.

Computations based on the model cannot solve this issue. Its a matter of scientific methodology. But computations can help to analyse this falsification of the model by putting under focus some variables which are more questionnable than others. If we are rather interested in model fitting, it could be usefull to know on what parameters to focus in order to make the minimal changes in the model when it does not agree with experimental data.

A solution to this issue can be given by considering the solutions of the model which are the closest to the experimental data. In discrete spaces as sign algebra cartesian products, the only interesting distance is the Hamming distance. Recall that the Hamming distance between two vectors (α 1 , . . . , α n ) and (β 1 , . . . , β n ) is the number of com- ponents such that α i = β i . So the Hamming distance between experimental data and solutions of E gives the minimal number of data which could be modified to fit the model. If variables represent parameters of the model, the Hamming distance gives the minimal number of adjustment of the model to fit the experimental data.

Of cause, it is not sufficient to have the minimum number of corrections. We are also interested in the closest vectors fitting the model and in the sets of minimal corrections. Let d h (v 1 , v 2 ) denote the Hamming distance between v 1 and v 2 two vectors of the same dimension on the sign algebra. If E is the set of solutions of E we define:

d h (v, E) = min w∈E (d h (v, w)) ( 13
)
The closest vector from v fitting the model is the Hamming projection of v onto E:

p h (v, E) = {w/ w ∈ E; d h (v, w) = d h (v, E)} (14)
and the set of minimal corrections is:

corr h (v, E) = {u/ v + u ∈ E & d h (v, v + u) = d h (v, E)} (15)

Computation of hard components

Computation of hard components is not a standard operation on polynomial functions on finite fields. However, it is possible using, Shannon decomposition as in standard operations. Let X = (x, X ) be a vector in ( / 3 ) n where the first component x corresponds to the first variable for the order used for QDD representation. The remaining components are components of the vector X . If

f (x, X ) = L α (0) f | 0 (X ) + L α (1) f | 1 (X ) + L α (-1) f | -1 (X ) (16)
is the Shannon decomposition of f then the computation of hard components is based on a few remarks: The representation of QDD by directed acyclic graphs result in large gains in memory compared to a raw representation as a tree derived from Shannon decomposition. However, all computations on QDD are based on this decomposition and thus potentialy have an exponential time complexity. The gain in memory is based on redondancy of many subtrees in the decomposition. This means that many computations are repeated if the algorithms are implemented by brute force application of Shannon decomposition. A fairly standard way to avoid repeated computations is to partially memorize them in a cache and this is used in the QDD package.

• (x, X ) is a solution of f (X) = 0 if and only if X is a solution of f | x (X ) = 0 • x is a hard solution of f (x, X) = 0 if
A technical problem with a cache is that it is generaly dedicated to cache one type of objects. The main cache of the QDD package is dedicated to QDD memorisation. We could build a cache for each type of object or some kind of universal cache. This would add a lot of complexity to the software. In our implementation we chosed to represent hard components as a polynomial function having the vector of hard components as unique root in the space of its variables. In this context, by variables of a polynomial, we mean variables appearing explicitely in the QDD representation.

Let X = (x, X ) and p is a polynomial function representing hard components among the variables of X . If x is a new hard component to be represented with value α, then a polynomial function representing the hard components of p and α is:

p = L α (x) × p (X ) + ∑ β =α L β × 1 (17)
With this in mind, it is easy to compute efficiently hard components of a polynomial function f considered as an equation. If

f (x, X ) = ∑ α∈K L α (x) f | α (X ) (18)
is the Shannon decomposition of f then the preceding remarks apply:

• if none of the f |α has roots then f has no roots and consequently no hard component.

• if only one of the f |α has roots, then x is a hard component with value the corresponding α. The hard components among the variables of X are the hard components of f |α.

• if several of the f |α have roots, then x is not a hard component and hard components of f are the hard components common to the f |α having roots and sharing the same value.

Hamming algorithms

In computing Hamming distance, projection and correction, we need to consider a more general case than the introductory one in 3.3. In general it is impossible to observe all the variables of a model. So the vector of observations V is indexed by a subset of the set of variables X = (X 1 , . . . , X n ). Let us introduce:

var(V ) = {X i /i is a component index of V } (19)
If the observations involve only a subset of variables, the Hamming distance from the set of solutions of E represented by a unique equation f (X) is redefined as:

d h (V, f ) = min{d h (V,W )/var(V ) = var(W ) & ∃Y f (V,Y ) = 0} (20)
It is the distance from V to the projection of the solutions of f (X) = 0 to var(V ).

Hamming correction definition has also to be made more precise:

corr h (V, f ) = {U/ var(V ) = var(U) & ∃Y f (V +U,Y ) = 0 & d h (V,V +U) = d h (V, f )} (21)
The Hamming projection is adapted in the same way;

pro j h (V ) = {U/ var(V ) = var(U) & ∃Y f (U,Y ) = 0 & d h (V,U) = d h (V, f )} (22)
.

In order to establish the recursive equations for the computation of the Hamming distance, Hamming projection and Hamming correction, we asssume that the variables are ordered and the components of V are ordered with the same order. By isolating the first variable we can write X = (x, X ) and V = [v,V ] t where exponent t denotes the transpose of a vector or a matrix. In the Shannon decomposition, we use the same order on variables:

f (x, X ) = L α (0) f | 0 (X ) + L α (1) f | 1 (X ) + L α (-1) f | -1 (X ) (23)
where at least two of the polynomial functions f | α are distinct so that x is the rootvar of f . The variable associated with the first component v of V is not necessarily equal to x. It might be equal, before or after x for the order on the variables. These three cases are the three cases of the recursive equation.

Identities on Hamming distance

We first establish simple identities on Hamming distance. Le E be a set of solutions of a system of qualitative equations E and W = [w,W ] t a vector in E. If V = [v,V ] t is a vector with var(v) = var(w), the following identity holds:

d h (V,W ) = d h (v, w) + d h (V ,W ) (24)
If true is coded by 1 and f alse by 0, d h (v, w) = (v = w). Taking the minimal value on W we get:

min W ∈E d h (V,W ) = min W ∈E (d h (v, w) + d h (V ,W )) (25) = min w∈E var(w) (d h (v, w) + min W ∈E w d h (V ,W )) (26)
where

E var(w) = {w/∃W s.t. (w,W ) ∈ E} (27) E w = {W /(w,W ) ∈ E} (28)
The first set is the projection of E onto the w component and the second one is the fiber over w. It is possible to get rid of the first set by setting min W ∈ / 0 d h (V ,W ) = ∞. With this convention we get:

min W ∈E d h (V,W ) = min w=0,1,-1 (d h (v, w) + min W ∈E w d h (V ,W )) (29
E w = {W / ∃Y f (w,W ,Y ) = 0} (30) = {W / ∃Y f | w (W ,Y ) = 0} (31)
and with d h (V, E) = d h (V, f ) when E is the set of solutions of the equation f (X) = 0 we can write:

min W ∈E w d h (V ,W ) = d h (V , f | w ) (32)
which leads to the recursive equation:

d h (V, f ) = min w=0,1,-1 (d h (v, w) + d h (V , f | w ) (33)
With the convention on empty sets we get the end cases for the recursion:

d h (V, 1) = d h (V, -1) = ∞ and d h (V, 0) = 0. Projection computation We assume that d h (V, f ) is computed. Let rh(V, f ) = {α/d h (V, f ) = d h (v, α)+d h (V , f | α )} and W an element of pro j h (V, f ). we get d h (V, f ) = d h (v, w) + d h (V ,W ). Now if Y is such that f (W,Y ) = 0, by Shan- non decomposition we have: f (W,Y ) = f | w (W ,Y ) = 0 (34) Moreover, d h (V , f | w ) = d h (V ,W ) since if W 1 were a vector such that there exist Y 1 with f | w (W 1 ,Y 1 ) = 0 and d h (V ,W 1 ) < d h (V ,W ), then W 1 = [w,W 1 ] t satisfies f (W 1 ,Y 1 ) = 0 and d h (V,W 1 ) < d h (V,W ),
showing that W is not in the Hamming projection pro j h (V, f ). We conclude that w ∈ rh(V, f ) and W ∈ pro j h (V , f | w ).

Reciproquely, if w ∈ rh(V, f ) and W ∈ pro j h (V , f | w ) then the vector W = [w,W ] t satisfies:

d h (V,W ) = d h (v, w) + d h (V ,W ) = d h (v, w) + d h (V , f | w ) = d h (V, f ) (35)
If we denote in the same way a polynomial having pro j h (V, f ) as solutions, we can write down a recursive equation on polynomial functions:

pro j h (V, f ) = ∑ w∈rh(V, f ) L w × pro j h (V , f | w ) + ∑ not w∈rh(V, f ) L w × 1 (36)

Correction computation

The computation of the Hamming correction derives directly from the previous equation. We look for vectors U such that V +U is in pro j h (V, f ). The preceding equation applied to U +V gives:

pro j h (V, f )(U +V ) = ∑ w∈rh(V, f ) L w (u + v) × pro j h (V , f | w )(U +V ) + ∑ not w∈rh(V, f ) L w (u + v) × 1 (37)
which is equal to zero if and only if u

+ v ∈ rh(V, f ) and pro j h (V , f | w )(U +V ) = 0 for w = u + v. This shows that U is in corr h (V , f | w-v
) and using the identity L w (u + v) = L w-v (u) we get:

corr h (V, f )(U) = ∑ w∈rh(V, f ) L w-v (u) × corr h (V , f | w-v )(U ) + ∑ not w∈rh(V, f ) L w-v (u) × 1 ( 38
)
again by using the same notation for a set and a polynomial equation representing this set.

Case 2: var(v) > rootvar( f )

Projection computation In this case, the polynomial function doesn't put any constraint on the components of vectors with index var(v). In other words, in the Shannon decomposition of f , with var(v) as first index,

f | 0 = f | 1 = f | -1 . In the Hamming distance expression d h (V, f ) = min w (d h (v, w) + d h (V , f | w )
, it is possible to have w = v and we get:

d h (V, f ) = d h (V , f | v ) = d h (V , f ) (39)
As an immediate consequence pro j h (V, f ) = pro j h (V , f ).

Correction computation

As the first component of vector V belongs to a vector in the Hamming projection, the correction on this component is u = 0. For the other part U of the vector U, we get from preceding relation:

d h (V +U, f ) = d h (V +U , f | v ) = d h (V +U , f ) (40)
showing that U is in corr h (V , f | v ). From these remarks we get for the polynomial function:

corr h (V, f )(U) = L 0 (u) × corr h (V , f ) + L 1 (u) × 1 + L -1 (u) × 1 (41)

Case 3: var(v) < rootvar( f )

This case includes the case V = / 0 which is not used directly but appears as an end case in recursion equations. We will first establish the recursive equations for V = / 0.

Distance computation

For computing the distance we start from the very definition (20):

d h (V, f ) = min Y {d h (V,Y )/ f (Y ) = 0} (42) = min Y {d h (V,Y )/ ∃y f (y,Y ) = 0} (43) = min Y {d h (V,Y )/ ∏ y=0,1,-1 f | y (Y ) = 0} (44) = d h (V, ∏ y=0,1,-1 f | y (Y )) (45) = d h (V , ∏ y=0,1,-1 f | y (Y )) (46) = min y=0,1,-1 d h (V , f | y ) (47)
since the set of solutions of ∏ y=0,1,-1 f | y (Y )) is the union of the sets of solution of each factor and

d h (V , f | y ) = ∞ whenever f | y has no root.
Projection computation From the preceding equations, if W satisfies d h (V, f ) = d h (V,W ) then there exists a Y and a y such that f | y (W,Y ) = 0 and the min is reached for y. So y ∈ rh(V, f ) and W is necessarely a vector of pro j

h (V, f | y ) . Reciproquely if y ∈ rh( f ,V ) and W ∈ pro j h (V, f | y ) then there exist a Y such that f | y (W,Y ) = 0. This vector satisfies the property: f (W, [y,Y ] t ) = f | y (W,Y ) = 0.
Since V has no variable corresponding to the index of rootvar( f ) which is the index of the component y, we get:

d h (V,W ) = d h (V, f | y ) = min y=0,1,-1 d h (V, f | y ) = d h (V, f ) (48)
since y ∈ rh(V, f ). so we get: pro j h (V, f ) = ∪ y∈rh(V, f ) pro j h (V, f | y ) and in polynomial function form:

pro j h (V, f ) = ∏ y∈rh(V, f ) pro j h (V, f | y ) (49)
Correction computation The computation of the correction derives directly from the computation of the projection. If V + U ∈ pro j h (V, f ) = ∪ y∈rh(V, f ) pro j h (V, f | y ) then there exists a y ∈ rh(V, f ) and V + U ∈ pro j h (V, f | y ). Consequently, U ∈ corr h (V, f | y ). The converse inclusion is obvious and we get:

corr h (V, f ) = ∏ y∈rh(V, f ) corr h (V, f | y ) (50)
as for the projection.

Case when V = / 0 Since this case appears only as an end case in recursive equations, let us consider the last step of a recursion ending with V = / 0. At this point we must have:

f = L 0 (x) f | 0 (X ) + L 1 (x) f | 1 (X ) + L --1(x) f |-1(X ) (51) V = [v] (52)
and var(V ) = {x}. On one hand we get from the recursive equation:

d h (V, f ) = min y=0,1,-1 (d h (v, y) + d h ( / 0, f | y )) (53) = min{d h ( / 0, f | v ), 1 + d h ( / 0, f | v 1 ), 1 + d h ( / 0, f | v 2 )} (54)
where v 1 , v 2 are the two other values in / 3 different from v.

On the other hand we have from our Hamming distance definition:

d h (V, f ) = min y=0,1,-1 {d h (v, y)/ ∃Y f | y (Y ) = 0} (55) = min f | y has roots = (d h (v, y)) (56) (57)
The only way to conciliate the two equations is to set:

d h ( / 0, f | α ) = 0 if f | α )has roots (58) d h ( / 0, f | α ) = ∞ either ( 59
)
Projection computation Since equation 36 must apply in the last step of the recursion, we are interested in the values which are in rh(V, f ). From equation (54) two cases appears:

• f v has roots Then v is the only element in rh(V, f ) and from pro j h (V, f ) = L v (x)pro j( / 0, f v ) + ∑ w =v L w (x)1 we must have pro j h (, f v ) = 0

• f v has no root Then the elements of rh(V, f ) are the w = v such that f w has roots. Since pro j h (V, f ) = L v (x)pro j( / 0, f v ) + ∑ w =v pro j( / 0, f w ) we must have pro j h ( / 0, f v ) = 1 and pro j h ( / 0, f w ) = 0 when f w has roots and 1 if not.

To summerize: pro j h ( / 0, f ) = 0 if f has roots (60) = 1 either (61)

Conclusion

We have presented original algorithms on polynomial functions on finite fields. These algorithms arised in the application of finite field to qualitative theory in order to provide efficient computation means for biological applications.

The algorithms have been implemented as a library written in C programming language. Python binding are also available providing easy access to these algorithms. The Python module, named Pyquali, has been sucessfully used in checking the coherence of data bases of gene interactions in E.Coli. A work under development uses the same kind of computation to infer a qualitative static model of gene interactions in Saccharomyces cerevisiae.

Figure 1 :

 1 Figure 1: From tree representation to direct acyclic graph for X 2 (Y + 1). The tree has 13 nodes while the DAG representing the same function has 5 nodes.

) 5 . 1

 51 Case 1: var(v) = rootvar( f ) Distance computation Applying the preceding identities on the set E = {W /∃Y f (W,Y ) = 0} we get:

  and only if there is only one value α for which f | α has roots.• a component x i of X is a hard component for f if and only if it is a hard component with the same value in every function f | α having roots..

In the implementation we use a dichotomic scheme which is equivalent but much more efficient

Correction computation Correction computation follows the same path. Using the same argument and starting from equation 38 and the expression of the distance (54) we get:

Random solution of a qualitative equation

In some applications we need to compute the proportion of solutions of a qualitative equation in a given reference space. Linked to this problem, is those of selecting a solution at random with a uniform distribution. If n is the number of solutions of a qualitative equation represented by a QDD f , we want to select a solution (x 1 , . . . , x k ) with probability 1 n . From Bayes's rule we get:

Let n α the number of solutions of f | α . Since (x 2 , . . . , x k ) is a solution of f | α drawn with uniform probability among all the solutions of f | α , we get:

Let p = n/3 k the proportion of solutions of the equation f in the space corresponding to the variables (X 1 , . . . , X k ) and p α = n α /3 (k-1) the proportion of solutions of f | α in the corresponding space. It is easy to establish the relation p = (p 0 + p 1 + p 2 )/3 from which we derive:

An algorithm may be builded on the preceding relations: first compute for each node of the QDD representing f , the proportion of solutions in the corresponding space. Then traverse this graph, down from the root, selecting each branch corresponding to the value α, with the probability given by (66). Unfortunately, such an algorithm would mean memorizing one float per node of the QDD, which may be quite hudge. To avoid this problem we propose to build a solution and compute the proportion of solutions at each node during the same traversal. The solution will be coded with the same trick we used for hard componants. This allows to use the general cache to memorize computations. Assume also that we have a function Random(π 0 , π 1 , π 2 ) which delivers a value α in the range {0, 1, 2} with the probabilities given as arguments.

Let p α be the proportion of solutions of f | α and sol α the random solution (coded as a QDD) drawn for f |α. If

is the Shannon decomposition of f , and ps = p 0 + p 1 + p 2 , we can draw α = Random(p 0 /ps, p 1 /ps, p 2 /ps) and then we get the recursive relations:

Stopping cases are obvious. Remark that a complete implementation of this algorithm must take into account the case where some variables associated with dimensions of the ambient space may not appear in the QDD f . This needs slight adaptations of the algorithm we will not describe here.