Saïd Amari

Isabel Demongodin

Jean Jacques Loiseau

Claude Martinez

Max-Plus Control Design for Temporal Constraints Meeting in Timed Event Graphs

Keywords: Discrete event systems (DES), feedback control, max-plus algebra, temporal constraints, timed event graph (TEG)

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

I. INTRODUCTION

Many manufacturing systems are subject to tight time constraints in their production process. Such constrained plants can be found in the semiconductor industry [START_REF] Kim | Scheduling analysis of time-constrained dual-armed cluster tools[END_REF], or in the automotive industry [START_REF] Martinez | Sizing of an industrial plant under tight time constraints using two complementary approaches: (max, plus) algebra and computer simulation[END_REF]. The aim of the presented work is to design a method to manage such a plant, so the strict time constraints are satisfied.

Production control and manufacturing plant problems are generally considered as Discrete Event Systems (DES) control problems. Petri Nets have been extensively used [START_REF] Cassandras | Introduction to Discrete Event Systems[END_REF], [START_REF] Murata | Petri nets: Properties, analysis and applications[END_REF] to deal with DES. Event graphs (a subclass of Petri nets) are often used, when the possible choices, for instance in the routing, are resolved with the help of the production control. We consider in the sequel a class of deterministic controlled processes modeled with Timed Event Graphs (TEG), that permit to take time constraints explicitely into account.

The problem tackled here could be formulated as a verification of some temporal conditions, see for instance [START_REF] Berthomieu | Modelling and verification of time dependant systems using time Petri nets[END_REF], [START_REF] Cofer | Control of event separation times in discrete event systems[END_REF], [START_REF] Ghezzi | A unified highlevel Petri net formalism for time critical systems[END_REF]. In the present contribution, the verification of time constraints is formulated in terms of a control problem, assuming that some inputs of the process can be controlled (it is generally the case). As is well-known from [START_REF] Baccelli | Synchronization and Linearity: An algebra for Discrete Event Systems[END_REF], TEG give rise to Max-Plus algebraic models, which are linear over the dioid max , and we make use of this linear framework. The plant behavior is modeled with Max-Plus equations, and the temporal constraints are represented with inequalities, also linear in the Max-Plus algebra. We propose a method for the synthesis of control law that permits to meet a given set of time constraints. The resulting control law itself is finally defined as a Max-Plus linear difference equation, involving a finite number of delays. This control law is causal, and can be implemented on-line, from the knowledge of the system state. Such an equation corresponds to a feedback that is also a TEG.

The control approach that we propose is quite different from that considered within the so-called supervisory control framework of timed discrete event systems ( [START_REF] Holloway | Asurvey of Petri net methods for controlled discrete event systems[END_REF], [START_REF] Moody | Feedback control of Petri nets based on place invariants[END_REF]). Time is explicitly taken into account in the proposed approach. TEGs and dioid framework have been used in the literature to treat a number of control problems. In [START_REF] Baccelli | Synchronization and Linearity: An algebra for Discrete Event Systems[END_REF] is calculated the latest control that permits to match a given target. The computation is based on residuation theory, and in terms of manufacturing processes, this control is interpreted in terms of a just-in-time production strategy. The residuation theory was further used by Cottenceau et al. [START_REF] Cottenceau | Synthesis of greatest linear feedback for timed event graphs in dioid[END_REF], to compute a feedback that delays as far as possible the input of the system without altering the input-output transfer relationship. It has also been used for time-varying systems in [START_REF] Lahaye | Just in time control of timevarying discrete event dynamic systems in (max; +) algebra[END_REF], and for uncertain systems with uncertain targets in [START_REF] Lhommeau | Interval analysis in dioid: Application to robust open-loop control for timed event graphs[END_REF]. Houssin et al. in [START_REF] Houssin | Just in time control of constrained (max; +)-linear systems[END_REF] propose a control strategy for the system to meet a set of general constraints, that includes time constraints. In the case that is considered here, the latest control (just-in-time) leads to a trivial behavior, if no more constraint is specified. The interpretation for a manufacturing system is that when the production is stopped, no part violates the time constraint, but this is not a pertinent solution to our problem. One seeks a control that is as permissive as possible, which means as early as possible. Unfortunately, there is no optimal control in that sense, since the usual method for the synthesis of control strategies for Max-Plus systems, that are based on residuation theory, does not apply in our case. An interesting formulation of our synthesis problem can be given in terms of controlled invariance. Katz [START_REF] Katz | Max-plus (A, B)-invariant spaces and control of timed discrete-event systems[END_REF] has shown that the satisfaction of linear inequality constraints defined in terms of a matrix E can be stated in terms of (A; B)-invariance of a semimodule included into the image of the Kleene star E 3 (see Section II-A for the definition of the Kleene star operation). The control design hence comes down to the computation of the maximal (A; B)-invariant semimodule included into the image of E 3 , as illustrated in Section VI of [START_REF] Katz | Max-plus (A, B)-invariant spaces and control of timed discrete-event systems[END_REF]. Unfortunately, this computation is difficult, and no general procedure is known for this computation. It may happen that this semimodule is not finitely generated, or that the (A; B)-invariant algorithm does not stabilize, because the dioid max is not artinian. We consider the case where not all components of the state are constrained, hence some entries of E 3 are not finite. As a consequence the condition proposed in [START_REF] Katz | Max-plus (A, B)-invariant spaces and control of timed discrete-event systems[END_REF] for the maximal (A; B)-invariant semimodule to be finitely generated is not satisfied, and we cannot use this formulation.

A first attempt to control TEG under strict temporal constraints has been presented in [START_REF] Amari | Control of linear min-plus systems under temporal constraints[END_REF]. This initial approach was developed in the Min-Plus algebra, under the assumption that all delays of the considered graph are integers. In the present contribution, this condition is not required, i.e., we consider a TEG with delays that may be real numbers.

The technical note is organized as follows. In Section II, some backgrounds on Max-Plus algebra and TEGs are recalled. The problem of finding a causal control law that satisfies critical time constraints is formulated in Section III. In Section IV, we first introduce a sufficient condition for the existence of a feedback under a single input and a single constraint. The method is constructive, and a suitable feedback is proposed, provided that the condition is satisfied. Then, the method is extended to the case of several temporal constraints. Finally, the general case, i.e., the characterization of a causal control law that guarantees several temporal constraints of a multi-input dynamic system modeled by a TEG is proposed in Section V. Section VI is devoted to an illustrative example. Concluding remarks are presented in Section VII.

II. BACKGROUNDS

A. Max-Plus Algebra

A monoid is a set, say D, endowed with an internal law, noted 8, which is associative and has a neutral element, denoted ". A semiring is a commutative monoid endowed with a second internal law, denoted , which is associative, distributive with respect to the first law 8, has a neutral element, denoted e, and admits " as absorbing element: 8a 2D, a " = " a = ". A dioid is a semiring with an idempotent internal law: 8a 2D, a 8 a = a. The dioid is said to be commutative if the second law is commutative. Max-Plus algebra is defined as ( [f01g; max; +). This semi-ring, denoted max , is a commutative dioid, the law 8 is the operation max with the neutral element " = 01, and the second law is the usual addition, with neutral element e =0 .Forn 2 and V , W 2 n2n max ,byV 8W we denote the matrix with components (V 8 W ) ij = max(V ij ;W ij ).Gi venp, q 2 , and matrices A 2 p2n max and B 2 n2q max , the symbol A B (or just A:B)

will denote the result of matrix multiplication defined by the formula

(A B) ij = n k=1 (A ik B kj ) = max k (A ik + B kj ):
The Kleene star of a square matrix M 2 n2n max , written M 3 is de- fined as M 3 = i2 M i , where M 0 equals the unit matrix, with entries equal to e on the diagonal, and " elsewhere. Let us recall that for x, v 2 n max the minimal solution of both the inequality, x M . x 8 v, and the equality, x = M x 8 v,i sx = M 3 v [START_REF] Baccelli | Synchronization and Linearity: An algebra for Discrete Event Systems[END_REF]. We denote by e i the vector of n max with its ith entry equal to e, and the other entries equal to ":

The following lemma and theorem will be useful to determine the existence of a solution to the problem of multivariable control tackled in this article. Lemma 1 concern a specific case of A x B x, see [START_REF] Allamigeon | The tropical double description method[END_REF] and [START_REF] Butkovič | An elimination method for finding all solutions of the system of linear equations over an extremal algebra[END_REF], [START_REF] Butkovič | Generators, extremals and bases of max cones[END_REF]. We provide it with its proof, that is simple, since it gives some hints for the proof of Theorem 1, which is more complicated, and for the sake of completion, because we use it in the proof of Lemma 2, in Section V.

Lemma 1: For any two row vectors u, v in 12n max , there exists a non trivial solution to inequality u x v x, if and only if there is k 2 n, such that u k v k .

Proof: The vector x = " is always a solution, called the trivial solution. A solution x is called non trivial if there exists an index k so that x k 6 = ". Under this condition, we can choose k so that also v x = v k + x k . Therefore we have u k + x k u x v x = v k + x k , which, since u k 6 = ", shows the necessity of the statement. Conversely, under the stated condition, one can choose x k = e, and x i = ", for i 6 = k, which defines a non trivial solution.

Theorem 1 gives explicitly the set of solutions to A x B x; details are given in [START_REF] Allamigeon | The tropical double description method[END_REF] and proof in [START_REF] Allamigeon | Computing the Extreme Points of Tropical Poly-Hedra[END_REF].

Theorem 1: Let C n max be defined as the set fx 2 n max j Ax Bxg, where A, B 2 q2n max (with q>0). Let G 0 ;G 1 ; ...;G q denote the sequence of finite subsets of n max defined as follows:

G 0 = e i ;i 2 n G `= fg 2 G `01 j A ` g B ` gg
[f(A ` h) g 8 (B ` g) h j g; h 2 G `01 ; A ` g B ` g; and A ` h>B ` hg for all `2 q, where A `and B `are the `-th rows of matrices A and B .

Then C is generated by the finite set Gq .

B. TEGs, Linear Max-Plus Models

An event graph is an ordinary Petri net where each place has exactly one upstream transition and one downstream transition. A TEG is obtained by associating delays to places or to transitions of a given event graph. In our case, delays are associated to places. Let t1 ; ...;tn be transitions having at least one upstream place, and let t u 1 ; ...;t u q be transitions having no upstream place and assumed to be controllable. If it exists, the place linking tj to ti , unique by hypothesis, is denoted by p ij , its corresponding delay is denoted ij and its marking, i. e. the number of tokens within this place, is denoted by mij .

Transition tj is controllable if there exists a path from transition t u s to transition t j . Such a path is an alternating sequence of transitions and places, of the form [t u s ;p k s ;t k ; ...;t k ;p jk ;t j].Bym we denote To represent the dynamic behavior of TEGs in Max-Plus algebra, a firing time us (k) of the kth occurrence is associated with each transition t u s and a firing time i (k) of the kth occurrence is associated with each transitions ti .

Example 1: The dynamic behavior of the TEG in Fig. 1 is characterized by the following equation:

3 (k) = max( 31 + u 1 (k 0 m 31 ); 32 + u 2 (k 0 m 32 ))
which, in Max-Plus algebra appears to be linear

3 (k)= 31 u 1 (k 0 m 31 ) 8 32 u 2 (k 0 m 32 )
The dynamic behavior of a TEG can be expressed by means of a linear equation in Max-Plus algebra as follows:

(k)= m0 [Am (k 0 m) 8 Bm u(k 0 m)] (1)
where the components of vector (k) are the firing times of the n transitions ti , the components of u(k) are the firing times of source transitions t u s . Matrix A m belongs to dioid n2n max , its entry A m equals to ij if there exists a place pij containing mij tokens, otherwise Am equals to ". Similarly, the entries of matrices Bm 2 n2q max correspond to the delays of places following source transitions. It is worth replacing (1) by the following explicit equation:

(k)= m>0 A 3 0 Am (k 0 m) 8 A 3 0 Bm u(k 0 m) (2)
where A 3 0 is the Kleene star of A 0 , previously defined (see [5] for more details). Entries of A0 are the delays associated to places without any token, therefore for a live event graph, entries of A 3 0 are such that

(A 3 0 ) ij < +1.
Analogously to the case of usual linear systems, the explicit (2) can be brought into a state space form. In order to obtain a state space model, one first expands all the places with a marking m>1 into m places with a marking equal to 1. Hence (m 0 1) intermediate transitions are added. Then the resulting extended state vector x(k) 2 N max is obtained, N = n +( m 0 1): The dynamic behavior of the expanded TEG is then described by an equation of the form, x(k)= Â0 x(k) 8 Â1 x(k 0 1) 8 B u(k), with Â0 , Â1 2 N2N max and B 2 N2q max which can be rewritten into the following explicit form, for k 1:

x(k)=A x(k 0 1) 8 B u(k) [START_REF] Amari | Sizing, Cycle Time and Plant Control Using Dioid Algebra[END_REF] where A = Â3 0 Â1, B = Â3 0 B, and initial conditions for `=1to N are x `(0) = ", the canonical initial conditions [START_REF] Baccelli | Synchronization and Linearity: An algebra for Discrete Event Systems[END_REF]. These formulations permit to point out that the behavior of a controlled TEG is fully determined by the input u(k), its state (3) and canonical initial conditions. Hence, the following formulation should be used:

x(k)=A ' x(k 0 ') 8 '01 k =0 A k B u(k 0 k 0 ) (4)
for each integer 1 ' k. Let us assume that the input u(k) is actually a control input, which can be arbitrarily assigned.

III. PROBLEM OF TEMPORAL CONSTRAINTS

Strict time constraints are frequent in industry. We took inspiration from the examples reported in [START_REF] Amari | Sizing, Cycle Time and Plant Control Using Dioid Algebra[END_REF] and [START_REF] Kim | Scheduling analysis of time-constrained dual-armed cluster tools[END_REF], that respectively concern rubber transformation for the automotive industry, and semiconductor manufacturing. One can for instance consider the example of a production process with a furnace for realizing a thermal treatment. The duration of any treatment in the furnace is fixed, or defined by a time interval. The aim is to control the system to meet this time constraint. A TEG already takes into account a delay on each place that corresponds to a minimal holding time. Nevertheless, in order to take into account the maximal duration, one has to express an additional constraint that the system should meet. The sojourn time of tokens in place p ij may be greater or equal to ij . In the presented work, a maximum sojourn time, noted max ij must also be respected. Hence, an interval of time [ij ; max ij ] is associated to place pij subject to a strict time constraint (see Fig. 2). This additional temporal constraint is expressed, for k mij, through the following inequality:

x i (k) max ij x j (k 0 m ij ):

(5)

IV. SINGLE CONTROL

In this section, TEGs with a single source transition t u are considered. Let us assume that the control law u(k) applied to this transition t u is given by a feedback u(k)=F x(k 0 1), for k 2 and u(1) = e. F is a row vector with N entries which are either Fr = " or F r 0, they correspond to the delay associated to control places that are to be determined to satisfy the constraints. Fr = " means the absence of control place between tr and t u . It also means that the control u does not depend on x r .

A. Single Constraint

Let us consider a TEG modelled by the linear Max-Plus (3) and subject to a single temporal constraint (5) on place p ij . Let us assume that transition t j is controllable, i.e., there exists a path from t u to t j where m is the cumulated marking along this path.

Theorem 2: Let be given the system (3), that defines the behavior of a TEG subject to a single constraint of the form (5) on place p ij . We assume that the cumulated marking along the path from the control input to place p ij is zero, m =0 , as well as the initial marking in the constrained place, m ij =0 . Then, the temporal constraint is satisfied applying a feedback of the form u(k)=F x(k01), if the following condition is satisfied:

Bi max ij Bj:

In addition, the following set of inequalities defines a suitable feedback, that guarantees that the temporal constraint is satisfied:

F r B j max ij A ir ; for r =1t oN: Fig. 3. TEG with temporal constraints on p and p . iterative procedure, as in Theorem 1. If the set is empty, no solution can be found to satisfy all Z constraints together. When a solution exists, a suitable feedback F should also be a solution of inequality [START_REF] Moody | Feedback control of Petri nets based on place invariants[END_REF].

Notice that if Z =1 , inequality [START_REF] Moody | Feedback control of Petri nets based on place invariants[END_REF] reduces to the condition of Lemma 2, which can be verified using Lemma 1.

VI. EXAMPLE

Let us consider the TEG of Fig. 3. This graph contains two source transitions modelling respectively u 1 (k) and u 2 (k) controls, q =2. The initial marking of place p 21 m 21 =0. There exists a path 1 from transition t u 1 to transition t1 and its initial marking is m =0.

One can check that max 21 B 11 B 21 , condition of Lemma 2, is satisfied. The initial marking of place p32 is m32 =0. There exists a path 2 from transition t u 2 to transition t2 and its initial marking is m =0.

One can check that max 32 B 2` B 3`; `2f1; 2g, i.e., condition of Lemma 2 is again satisfied. Furthermore, choosing F `r = " 8`, r 2 N except F 13 = F 14 = F 15 = F 24 =0, results in a suitable feedback Fig. 5. Resulting controlled TEG: temporal constraints on p and p are satisfied by places added between transitions t and t , t and t , and transitions t and t . that satisfy max 21 2 `=1 (B 1` F `r ) A2r 8 2`=1 (B 2` F `r ), for r =1to 5, and max 32 2 `=1 (B 2` F `r ) A 3r 8 2 `=1 (B 3` F `r ),for r =1to 5. Finally, according to Theorem 4, the control law which guarantees the satisfaction of both temporal constraints is given by u(k)= x3(k 0 1) 8 x4(k 0 1) 8 x5(k 0 1)

x 4 (k 0 1)

, that can be reduced to u(k)= x 3 (k 0 1) 8 x 4 (k 0 1) x4(k 0 1)

. This feedback can be interpreted by three control places connected to the TEG to guarantee the respect of the temporal constraints. The controlled graph is given in Fig. 5.

VII. CONCLUSION

We have proposed sufficient conditions for the existence of a control for TEGs subject to strict temporal constraints. The results include the case multivariable control. This method is illustrated on a simple example from manufacturing. The method is explicitely constructive, provided that the conditions hold, leading to a simple control design.

The presented approach in the Max-Plus algebra resembles to the one in the Min-Plus algebra already developed in [START_REF] Amari | Control of linear min-plus systems under temporal constraints[END_REF]. We point out that the conditions obtained in Max-Plus algebra are simpler than that obtained in the Min-Plus case. Just one inequality has to be checked in Theorem 2, which is a very simple test to the existence of a solution. In addition, the Max-Plus control law is easier to compute, since there are no need to compute matrix exponent as in the Min-Plus case. Furthermore, depending on the plant to be studied, it can be a great advantage to deal with the Max-Plus approach, especially if the delays involved in the model lay within a wide range of values: this could lead to a very large state vector in the Min-Plus case.

One can note that both Min-Plus and Max-Plus approaches lead to a supervisor that, if the conditions are fulfilled, guarantees that the time constraints are satisfied. For a system with an empty path (m =0) from the control transition t u to the input transition t j of the con- strained place, a trivial control law u(k)=+1, for k 2 guarantees that the time constraints are satisfied. In a production plant context, this trivial control is unacceptable, and one searches for a control law that does not unnecessarily slow down the manufacturing plant throughput. The proposed approach gives a minimal supervisor, which in general is not infimal, since an infimal supervisor does not always exist. This contrasts the work of [START_REF] Cottenceau | Synthesis of greatest linear feedback for timed event graphs in dioid[END_REF], [START_REF] Houssin | Just in time control of constrained (max; +)-linear systems[END_REF], [START_REF] Lahaye | Just in time control of timevarying discrete event dynamic systems in (max; +) algebra[END_REF]- [START_REF] Maia | Optimal closed-loop control for timed event graphs in dioid[END_REF], where supremal supervisors are proposed, to solve various control problems, that do refer to Just-In-Time optimality criteria.

It was noticed in the introduction that an alternative formulation for this family of problems is in terms of (A; B)-invariance. This fact was pointed out by Katz [START_REF] Katz | Max-plus (A, B)-invariant spaces and control of timed discrete-event systems[END_REF], and generated recent work, see for instance [START_REF] Di Loreto | Duality between invariant spaces for max-plus linear discrete event systems[END_REF]. The solution in this context would come from the computation of a maximal (A; B)-invariant semimodule, which is an open problem in general (see [START_REF] Katz | Max-plus (A, B)-invariant spaces and control of timed discrete-event systems[END_REF]). The present material may help to construct new examples toward a solution to this difficult and important problem.

Fig. 1 .

 1 Fig. 1. TEG example.

Fig. 2 .

 2 Fig. 2. Temporal constraint.

Fig. 4 .

 4 Fig. 4. TEG of Fig. 3 with additional transition t .

  Two additional temporal constraints f max 21 ; max 32 g are assigned to places p21 and p32 of this graph, and they are expressed respectively by the inequalities x 2 (k) 1 x 1 (k), and x (k) 1 x (k). The problem consists in calculating a control vector u(k)=( u1(k)u2(k)) t , which satisfies both constraints. The graph of Fig.3has been transformed into the graph of Fig. 4. To do so, place p 12 containing 2 tokens has been split into two places marked to 1 and the intermediate transition t5 has been added.The state equation associated with this new TEG is where the components of x(k) are the firing times of transitions t 1 , t 2 , t 3 , t 4 and t 5 , and the vector u(k) is the control law. The meeting of constraints

ACKNOWLEDGMENT

The authors wish to thank the associate editor and the anonymous reviewers, for greatly helping to improve this technical note, that has evolved a lot from the initial version.

Proof: First remark that, by assumption, there exists a path from the control transition t u to transition tj, then Bj 6 = ". Hence, a finite feedback F will always exist under the conditions of Theorem 2. According to (4), the ith component of the vector x(k) is given by the following explicit expression: [START_REF] Berthomieu | Modelling and verification of time dependant systems using time Petri nets[END_REF] for every integer ' 1. Then, considering only the contribution of the control input, e.g., the relation between transitions t u and t j , we obtain that the jth state vector component x j (k) satisfies [START_REF] Butkovič | An elimination method for finding all solutions of the system of linear equations over an extremal algebra[END_REF] (A m B) j u(k 0 m ) x j (k): [START_REF] Butkovič | An elimination method for finding all solutions of the system of linear equations over an extremal algebra[END_REF] Taking [START_REF] Berthomieu | Modelling and verification of time dependant systems using time Petri nets[END_REF] into account, we obtain that the constraint ( 5) is satisfied if the two following inequalities hold:

Further, taking [START_REF] Butkovič | An elimination method for finding all solutions of the system of linear equations over an extremal algebra[END_REF] into account, we obtain that inequalities ( 8) and ( 9) are satisfied if the following inequalities hold:

As the function u(k) is non-decreasing, the inequality ( 14) is satisfied if the inequality k 0 m + m ij is satisfied for k 0 =0to ' 0 1.B y assumption m = mij =0this inequality is always true. Again with m = m ij =0 ,wehave' =1and hence condition [START_REF] Di Loreto | Duality between invariant spaces for max-plus linear discrete event systems[END_REF] reduces to F r B j max ij A ir , for r =1to N , condition (13) reduces to Bi max ij Bj.

Remark 1: The hypotheses that m = m ij =0 , are met in general for production plants. It means that at initial state, there is no product in process in the plant, which is not restrictive in practice.

Remark 2: A suitable feedback F is found choosing Fr = " for A ir = " and F r =max(0;A ir 0 B j 0 max ij ) for A ir 6 = ".

B. Multiple Constraints

We consider now the case of a TEG having one source transition, which is controllable, with Z constrained places. These places are denoted p z , for z =1to Z. For each constrained place p z , let m z , z and max z respectively denote the initial marking, the minimal and maximal delays. Further, let tz and t z respectively denote the input and output transitions of place p z , let x z (k) and x z (k) denote the corresponding firing dates, and let m denotes the cumulated marking along path z going from source transition t u to tz. The temporal constraints are now expressed by the inequalities

for z =1to Z. We denote by u z (k) the control law calculated as in the previous section to satisfy the zth temporal constraint. The following theorem defines a causal feedback, which ensures that all Z temporal constraints are satisfied.

Theorem 3: Let be given the system (3), subject to Z time constraints of the form [START_REF] Holloway | Asurvey of Petri net methods for controlled discrete event systems[END_REF], where mz = m =0, for z =1to Z. Then, there exists a causal control law of the form u(k)= Z z=1 u z (k), with u z (k)=Fz x(k01), that ensures the satisfaction of the Z temporal constraints, if the following condition is satisfied: to Z, it is clear that u(k)= Z `=1 u `(k) fulfills all the Z temporal constraints, which ends the proof.

V. M ULTIVARIABLE CONTROL

Let us consider in this section a TEG with q source transitions, q 1. The behavior of such a graph is again represented by a linear Max-Plus system (3), the control law being a vector u of q components. Let us first suppose that a single place p ij is subject to a temporal constraint of the form [START_REF] Baccelli | Synchronization and Linearity: An algebra for Discrete Event Systems[END_REF]. The problem is to design a control law u(k) 2 q max , with q 1, to satisfy the constraint [START_REF] Baccelli | Synchronization and Linearity: An algebra for Discrete Event Systems[END_REF]. The components of u(k) are denoted by u `(k) for `=1to q. We denote by m the cumulated markings along path `from a controllable transition t u `to transition t j .

Lemma 2: Let a TEG with q source transitions (q 1) be subject to a single temporal constraint of the form (5) on place pij. Then, there exits a feedback of the form u(k)=F x(k 0 1) to guarantee the satisfaction of the temporal constraint (5) if there exists an index `such that max ij B j` B i`, with m =0.

A suitable feedback is obtained with

for r =1to N.

Proof: If all B j`= ", no solution can be found. This means that no path exists from any controllable transition t u `to transition t j .I f there exist q controllable transitions, then the jth state vector component xj(k) satisfies the inequality q l=1 (A m B) j` u `(k 0 m ) xj(k): [START_REF] Katz | Max-plus (A, B)-invariant spaces and control of timed discrete-event systems[END_REF] According to (4), the ith component of vector x(k) is given by the following explicit expression:

for every integer ' 1. Taking ( 17) into account, we obtain that the constraint ( 5) is satisfied if the two following inequalities hold:

Further, taking ( 16) into account, inequalities ( 18) and ( 19) are satisfied if the following inequalities hold: `=1 B i` F `r x r (k) q `=1 max ij B j` F `r xr(k) for r =1to N, hence if q `=1 B i` F `r q `=1 max ij B jl F `r for r =1to N. As stated by Lemma 1, if 9`2 q j max ij B j` B i`, then there exists a solution vector F r for each column r.

Let us now consider a TEG with Z temporal constraints. If there exists a path from a controllable transition t u s to transition t z for each Z constraints, with m = mz =0, then one can synthesize a control law that satisfies all Z temporal constraints using Theorem 4.

Theorem 4: Let a TEG with q source transitions (q 1) be the then the Z constraints are met using the causal control law u(k)= Z z=1 u z (k), where u z (k) is the control law given by Lemma 2, that guarantees that the zth constraint is satisfied.

Proof: A control law u z (k) satisfies the zth constraint if and only if Lemma 2 holds. Thus, in order to satisfy condition of Lemma 2 for each z =1to Z constraints, one has to determine a feedback F that satisfy q `=1 B z ` F `r q `=1 max z B z l F `r for r =1 to N, and for all s =1to Z. These Z conditions are summarized in inequality [START_REF] Moody | Feedback control of Petri nets based on place invariants[END_REF]. The set of solutions to such an inequality is given by an


Saad Sheikh

Rolf Backofen

Yann Ponty

Impact Of The Energy Model On The Complexity Of RNA Folding With Pseudoknots

Keywords: RNA folding, General pseudoknots, Hardness, Inapproximability

Predicting the folding of an RNA sequence, while allowing general pseudoknots (PK), consists in finding a minimal free-energy matching of its n positions. Assuming independently contributing basepairs, the problem can be solved in Θ(n 3 )-time using a variant of the maximal weighted matching. By contrast, the problem was previously proven NP-Hard in the more realistic nearest-neighbor energy model. In this work, we consider an intermediate model, called the stackingpairs energy model. We extend a result by Lyngsø, showing that RNA folding with PK is NP-Hard within a large class of parametrization for the model. We also show the approximability of the problem, by giving a practical Θ(n 3 ) algorithm that achieves at least a 5-approximation for any parametrization of the stacking model. This contrasts nicely with the nearest-neighbor version of the problem, which we prove cannot be approximated within any positive ratio, unless P = N P .

Introduction

Ribonucleic Acid (RNA) is one of the key pieces to the puzzle of molecular biology. It plays a very large number of roles, not only by coding for proteins, but also through catalytic and regulatory functions. To play such roles, RNA folds into an intricate structure which is stabilized by the pairing, mediated by hydrogen bonds, of some of its positions. The conformations that arise from this folding process are instrumental to the function of an RNA. Consequently, the process of RNA folding has been extensively studied by molecular biology and biochemistry, and its in silico prediction has given rise to a wealth of computational approaches. Early work on the subject have focused on the secondary structure, a restriction of all admissible base-pairs that forbids crossing-interactions. Under Fig. 1. Canonical base-pairs of the Oceanobacillus iheyensis group II intron, derived from a 3D model (PDB id: 3IGI) using RNAView [START_REF] Yang | Tools for the automatic identification and classification of rna base pairs[END_REF].

the assumption of reasonable, additive, energy models, such a restriction implies an optimal-substructure property for computing the most stable conformation, i.e. the one having minimal free-energy. Polynomial-time algorithms, based on dynamic-programming (DP), have consequently been proposed for predicting the minimal free-energy secondary structure conformation of an RNA from its sequence. Unfortunately, the assumption of a non-crossing conformation may impede the quality of the prediction, since functionally essential crossing-interactions are found in many families of non-protein-coding RNAs (ncRNAs), as illustrated by Figure 1. For instance, pseudoknots (PK), can be found within the RFAM consensus [START_REF] Griffiths-Jones | Rfam: an RNA family database[END_REF] of at least 70 functional families of ncRNAs and are conserved throughout the evolution.

Taking general pseudoknots into account is known to turns RNA minimal-free energy prediction into a rather challenging problem. A pioneering work by Cary, Tabaska et al [START_REF] Cary | Graph-theoretic approach to RNA modeling using comparative data[END_REF][START_REF] Tabaska | An RNA folding method capable of identifying pseudoknots and base triples[END_REF] considered a simple additive model, associating independent energy contributions to each putative base-pair, and used a O(n 3 ) maximal-weighted matching algorithm to extract a minimal free-energy folding. Unfortunately, this energy model is regarded as unrealistic because of its incapacity to capture the interaction of consecutive -stacking -base-pairs, which constitute a primary stabilizing force in RNA folding. Such energy contributions are captured by the nearest neighbor energy model, in which the contribution of each base-pair depends on the base-pairing status and partners of its consecutive positions. The hardness of RNA folding assuming a nearestneighbor energy model was independently established by Lyngsø and Pedersen [START_REF] Lyngsø | RNA pseudoknot prediction in energy-based models[END_REF], and Akutsu [START_REF] Akutsu | Dynamic programming algorithms for RNA secondary structure prediction with pseudoknots[END_REF]. Subsequent efforts have therefore focused on providing either parameterized complexity algorithms [START_REF] Liu | RNA folding including pseudoknots: A new parameterized algorithm and improved upper bound[END_REF][START_REF] Zhao | Rapid ab initio prediction of RNA pseudoknots via graph tree decomposition[END_REF], heuristics [START_REF] Bindewald | Cylofold: secondary structure prediction including pseudoknots[END_REF][START_REF] Theis | Prediction of RNA secondary structure including kissing hairpin motifs[END_REF] or exact DP schemes for tractable subsets of pseudoknots [START_REF] Rivas | A dynamic programming algorithm for RNA structure prediction including pseudoknots[END_REF][START_REF] Reidys | Topology and prediction of RNA pseudoknots[END_REF].

It is frequent that the complexity of solving any problem in computational biology optimally is tied to the chosen model (e.g. [START_REF] Ashley | On approximating an implicit cover problem in biology[END_REF]). However, despite a significant amount of research focusing on predicting pseudoknots, the impact of a specific instantiation of the energy model on the computational complexity of RNA folding with pseudoknots has only been partly unexplored. In this extended abstract, we further study the influence of the energy model on the complexity and approximability of RNA folding with unconstrained pseudoknots. In addition to the base-pair and nearest-neighbor models, we consider the stacking base-pairs energy model, which captures the dependency between consecutive base-pairs. The computational complexity of RNA folding under this energy model was first studied by Ieong et al [START_REF] Ieong | Predicting RNA secondary structures with arbitrary pseudoknots by maximizing the number of stacking pairs[END_REF]. They were able to show the NPhardness of maximizing the number of stacking-pair among the set of planar secondary structures, a restriction of general pseudoknots. This restriction was lifted by Lyngsø [START_REF] Lyngsø | Complexity of pseudoknot prediction in simple models[END_REF], who established the hardness of maximizing the number of base-pairs, allowing general types of pseudoknots. Approximation algorithms we also sought, leading to the current best 8/3-approximation O(n 10 )-time algorithm reported by Jiang [START_REF] Jiang | Approximation algorithms for predicting RNA secondary structures with arbitrary pseudoknots[END_REF]. However all of these works consider a purely combinatorial model, maximizing the number of base-stacking, while the contribution of stacking pairs to the free-energy may vary significantly. It is therefore a natural question to ask to what extent the hardness of folding with pseudoknots is affected by perturbations of the energy model. More generally, understanding what makes the problem hard, and just how hard, could be instrumental to the development of future algorithms, achieving better tradeoffs between sensibility and complexity.

This extended abstract is organized as follows. First we formally define in Section 2 our main problem, along with the different energy models considered. We discuss the NP-hardness of the stacking base-pairs in Section 3, and present an approximation in Section 4. In Section 5, we show the inapproximability of RNA folding with pseudoknots under the nearest-neighbor energy model. Finally Section 6 summarizes the contributions and describes futures lines of research.

Problem statement and free-energy models

Let ω ∈ {A, C, G, U} * be an RNA sequence, and m be a partial matching of the positions in ω, i.e. a set of non overlapping pairs of positions in ω. An energy model is a real-valued function E w that associates a freeenergy to ω by summing over the contributions of local motifs in the matching. The precise definition of local motifs will depend on the exact energy model. A low free-energy indicates a stable folding. Furthermore, any freeenergy contribution is usually determined up to an additive constant. Therefore one can assume that the contribution to the free-energy of any local motif is negative, with the exception of extremely unfrequent motifs which will be forbidden and assigned +∞ contributions. Let us then rephrase the problem of as an optimization (minimization) problem.

RNA-PK-Fold(E) problem

Input: An RNA sequence w. Ouput: A partial matching m over w, i.e. a set of pairwise disjoint pairs of positions in [1, |w|], which minimizes E w (m).

The three reference energy models are usually considered:

-Base-pairs model B [START_REF] Nussinov | Fast algorithm for predicting the secondary structure of single-stranded RNA[END_REF][START_REF] Cary | Graph-theoretic approach to RNA modeling using comparative data[END_REF][START_REF] Tabaska | An RNA folding method capable of identifying pseudoknots and base triples[END_REF]: Here, local motifs are simply individual base pairs, independently contributing to the free-energy:

B w (m) = (i,j)∈m ∆ B (w i , w j )
where

∆ B : {A, C, G, U} 2 → R -∪ {+∞}.
-Stacking pairs model S [START_REF] Lyngsø | Complexity of pseudoknot prediction in simple models[END_REF][START_REF] Bon | Prédiction de structures secondaires d'ARN avec pseudo-noeuds[END_REF]: This model only considers consecutively nested pairs as motifs, and disregards isolated pairs:

S w (m) = (i,j),(i+1,j-1)∈m ∆ S (w i , w j , w i+1 , w j-1 )
where

∆ S : {A, C, G, U} 4 → R -∪ {+∞}.
-Nearest-neighbors model N [START_REF] Lyngsø | RNA pseudoknot prediction in energy-based models[END_REF][START_REF] Rivas | A dynamic programming algorithm for RNA structure prediction including pseudoknots[END_REF]: This motif definition is even more expressive, allowing different contributions for each base-pair, depending on its bases, the base-pairing status of its consecutive neighbors and own their own partners:

N w (m) = (i,j)∈m i<j ∆ N (w i , w j , w i+1 , w j-1 , w m i+1 , w m j-1 )
where

∆ N is any function {A, C, G, U} 4 × {A, C, G, U, ∅} 2 → R -∪ {+∞}, m i denotes the partner of a position i in m (or ∅ if i is unpaired, while w ∅ ≡ ∅ by convention).
These three models induce different decompositions into motifs for any given structure, as illustrated by Figure 2.

NP-hardness of RNA-PK-Fold(S) in any non-degenerate stacking energy model

Consider the set of canonical base-pairs (A, U), (G, C) and (G, U).

A combinatorial stacking model S * specializes the stacking pairs model by assigning a -1.0 kcal.mol -1 contribution to any canonical stacking pair, and +∞ to others. It was showed by Lyngsø [13] that the RNA-PK-Fold(S * ) problem is NP-complete, using a reduction from the BIN-PACKING problem.

Here we complement this result by showing its robustness, i.e. the NP-hardness of the problem under a wide class of stacking energy model.

Theorem 1 Let S be a stacking energy model that allows (G, C) pairs, and forbids (A, C) and (A, G) pairs. Then RNA-PK-Fold(S) is NP-hard.

Proof. In order to prove the hardness of RNA-PK-Fold(S), let us remind the statement of the 3-PARTITION problem:

3-PARTITION problem Input: A multiset of integral values X = {x i } n i=1 of cardinality n = 3m, such that n i=1 x i = m • K for some K ∈ N, and ⌊K/4⌋ < x < ⌊K/2⌋, ∀x ∈ X. Ouput: True if there exists a partition of X into m triplets ((x a j , x b j , x c j )) m j=1 such that x a j + x b j + x c j = K, ∀j ∈ [1, m],
and False otherwise. From Garey and Johnson [START_REF] Garey | Complexity results for multiprocessor scheduling under resource constraints[END_REF], it is known that 3-PARTITION is strongly NP-complete, i.e. not only is the problem NP-hard, but it remains hard even when the elements of X are upper-bounded by some polynomial function P (n).

A B C D X =            2 4 2 3 2 1 x1 x2 x3 x4 x5 x6            (K = 7) C C A C C C C A C C A C C C A C C A C A G G G G G G G A G G G G G G G wX = G G G G G G G A C C A C A A C C A C C C C A C C A C C C G G G G G G G C x1 C x2 C x3 C x4 C x5 C x6 G K G K
Lemma 2 Let X be a 3-PARTITION instance whose values are bounded by P (n), and w X be an RNA sequence such that:

w X = C x 1 AC x 2 AC x 3 A • • • AC xn AG K AG K A • • • AG K m times
.

There exists a 3-partition of X into equally summing triplets if and only if there exists a solution to RNA-PK-Fold(S) over w X having energy

k = δ • (K -3) • m kcal.mol -1 with δ := ∆ S (C, G, C, G).

Let us summarize the argument:

-A matching has minimal free-energy iff any block C x is entirely paired to some contiguous substring of a single G K block. -A matching has minimal free-energy iff every position in every G K block is connected. -Any optimal matching thus gives us a mapping between the C x and G K blocks, which can be transformed in polynomial-time into a solution to the 3-Partition problem.

Proof. X is 3-partitionable ⇒ ∃m * such that S w X (m * ) = δ•(K-3)•m:
If X is 3-partitionable, then there exists m disjoint triplets ((x a j , x b j , x c j )) m j=1 whose sum is identically K. If follows that the C x blocks in w X can be partitioned into triplets (C xa j , C x b j , C xc j ) that can form a three-tier perfect helix with the j-th block G K ≡ G xa j .G x b j .G xc j . By creating x a j (resp.

x b j and x c j ) nested base-pairs between the block C xa j (C x b j and C xc j ) and the beginning (resp. middle and ending) of the j-th block G K , one obtains exactly (x a j ) -1 + (x b j ) -1 + (x c j -1) = K -3 stacking pairs. Repeating the operation for each triplet j yields a valid conformation with (K -3) • m stacking (G, C)/(G, C) pairs, and the implication follows.

∃m * such that S w X (m * ) = δ • (K -3) • m ⇒ X is 3-partitionable:
Let us remark that the absence of U implies that the only admissible base-pairs are C/G or G/C, arising from interactions between C x and G K blocks respectively.

First, let us show that each G K block contributes to at most K -3 stacking pairs, and that this upper-bound cannot be reached unless G K creates K base-pairs with exactly 3 distinct C x blocks. Indeed, it is easily seen that any G K block, connected to b blocks (C x d 1 , . . . , C x d b ) by at least one base-pair, for a total number P of base-pairs, creates at most Pb stacking pairs. This bound is reached when G K is split into b portions, each forming a perfect helix with the corresponding C x d i block. Noting that x i < ⌊K/2⌋ is equivalent to x i ≤ ⌊K/2⌋ -1, one has that any connection of G K with b blocks can therefore create at most min(K, b • (⌊K/2⌋ -1)) base-pairs. It follows that, for b = 1 and b = 2, the maximum number of stackings involving G K is bounded by ⌊K/2⌋ -2 and 2⌊K/2⌋ -4 ≤ K -4 respectively. For b ≥ 3, the number of base-pairs is potentially no longer limited by the lack of occurrences of C, but by the K occurrences of G in G K . It follows that, when b ≥ 3, the maximum number of stacking pairs is K -3, and is reached for b = 3 when every position in G K is paired.

Then let us assume the existence of a matching with (K -3) • m stacking pairs. Since K -3 is the upper-bound on the number of stacking pairs supported by a given G K block, then each of the m G K blocks must achieve this upper-bound. It follows that each G K block must create a total of exactly K base-pairs with a triplet of blocks (C x a , C x b , C x c ). A direct corollary is that every G and C in w X must be paired.

We have now established that, within any matching having m • (K -3) stacking pairs, each G K block creates exactly K base-pairs with a triplet of blocks (C

x a i , C x b i , C x c i ).
To conclude on the implication, we need to show that each C x block interacts with a single G K block, i.e. that the (C x a i , C x b i , C x c i ) triplets are mutually disjoint. Indeed, if a block C x i is found in two distinct triplets, then there exists a block C x j that is not within any triplet (remember that there are 3K blocks C x and K triplets). It follows that at most m•K-x j base-pairs exist within this matching, which contradicts K base-pairs for every G K block. Consequently, no C x i can be present in two distinct triplets, and the triplets are therefore disjoint. Therefore, the interacting blocks found in a matching having energy δ • (K -3) • m induce a partition of the {C x i } 3m i=1 blocks into m triplets. Furthermore, each (C x a , C x b , C x c ) triplet must give rise to K base-pairs, and therefore x a + x b + x c ≥ K. Since the triplets are disjoint and partition a set of a total m • K occurrences of C, then any excess of C within a triplet implies a lack of C within another, so one has x a + x b + x c = K. We conclude that any matching of w X having energy δ • (K -3) • m induces the existence of a partition {C x i } 3m i=1 blocks into disjoint triplets (C xa , C x b , C xc ) such that x a + x b + x c = K which, in turn, implies the existence of a 3-partition for X.

⊓ ⊔

It follows from Lemma 2 that any algorithm for RNA-PK-Fold(S) gives an algorithm for the 3-PARTITION problem. Furthermore the length of w X exactly equals n i=1

x i +K•m+2m-1 = 2K•m+2m-1 ∈ O(n 2 •P (n))
where P (n) is the polynomial upper bound on the value of each x i . Therefore any polynomial algorithm for RNA-PK-Fold(S) gives a polynomial algorithm for the 3-PARTITION problem. Since 3-PARTITION is NPhard, then so is RNA-PK-Fold(S) and Theorem ?? follows.

⊓ ⊔

Approximability of RNA-PK-Fold(S) in the stacking model

Since objective functions are usually derived experimentally or statistically, it is a natural question to ask whether hard problems can be efficiently approximated. Previous works on the subjetc only considered a combinatorial version of the problem, and the current best algorithm [START_REF] Jiang | Approximation algorithms for predicting RNA secondary structures with arbitrary pseudoknots[END_REF] produces a matching whose number of stacking pairs is guaranteed to be at least 3/8 • OP T , where OP T is the maximal number of stacking pairs in any matching. Unfortunately, this result does not hold for arbitraryvalued stacking energy models, as the free-energy of valid stacking pairs may greatly vary. For instance, the latest version of the Turner model reports a factor ∼ 3.6 discrepancy between stacking canonical pairs, bringing the guaranteed approximation ratio down to 1/10. By contrast, we show that RNA-PK-Fold(S) can be approximated in polynomial time up to a factor at least 1/5, for any stacking model S.

Input : An RNA sequence w Output: A matching m of non-overlapping pairs of positions

G = (V, E) ← ([1, |w| -1], ∅); M ← ∅; foreach u, v ∈ V do
if wu base pairs with wv+1 and wu+1 base pairs with wv then // Label each edge with its weight/energy E ← E ∪ (u, v, -∆S (wu, wv+1, wu+1, wv));

end end m ′ ← MaxWeightedMatching(G); foreach (u, v) ∈ m ′ sorted by increasing value ∆S (wu, wv+1, wu+1, wv) do if ∀(u ′ , v ′ ) ∈ m, {u ′ , v ′ } ∩ {u, v + 1, u + 1, v} = ∅ or (u ′ , v ′ ) ∈ {(u, v + 1), (u + 1, v)} then m ← m ∪ {(u, v + 1), (u + 1, v)}; end end return m
Algorithm 1: A 5-approximation for any stacking energy model.

Theorem 3 In any stacking energy model, RNA-PK-Fold(S) ∈ AP X, and can be approximated in polynomial time within a factor at least 1/5.

Proof. To prove the approximability of RNA-PK-Fold(S), let us consider Algorithm 1. This algorithm contracts consecutive positions in the RNA sequence as vertices, and adds an edge, weighted according to the energy function, between any pair of compatible positions. Computing a maximal weighted matching on this graph gives a set of stacking-pair which is not necessarily a valid matching, since distinct pairs of stacking pairs may induce more than a single partners for a given position. Therefore the algorithm considers the returned stacking pairs in decreasing order, and only retains the stacking pairs that do not conflict with the current selection of stacking-pairs. Now let m * be the optimal matching for the given RNA string w, m ′ be the maximal matching over G, and m be the matching finally returned by the algorithm. Let us remark that m ′ induces a set of matched pairs over w that does not strictly constitutes a matching, as some position may be matched twice. Nevertheless let us write S w (m ′ ) as a shorthand for the total energy of m ′ , obtained by summing over the stacking pairs induced by m ′ . Any matching can be decomposed as a set of stacking pairs (leaving a set of isolated, non-contributive, base-pairs), hence one has S w (m ′ ) ≤ S w (m * ) ≤ 0. Any edge (i, j) in m ′ may conflict with at most 4 other, adjacent, stacking-pairs. Furthermore, the algorithm con-siders the edges in m ′ by decreasing contribution, so the stacking pairs induced by any edge (i, j) in m ′ may only conflict with four stacking pairs having (negative) contribution of smaller absolute value. Discarding these competitors guarantees that at least 1 5 of the total energy of m ′ is retained in m, i.e. S w (m) ≤ 1 5 S w (m ′ ) ≤ 0, and one therefore concludes that Sw(m) Sw(m * ) ≥ 1 5 .

⊓ ⊔

Remark that the actual approximation ratio achieved by Algorithm 1 might be better than 1/5, even in the worst-case scenario. However, this crude upper-bound already establishes the approximability of the problem, nicely contrasting with our upcoming inaproximability result for the nearest-neighbor version of the problem.

5 Inapproximability of RNA-PK-Fold(N ) in the nearest-neighbor energy model

The stacking model, considered in the above sections, makes the prediction of RNA structure NP-Hard, yet remains approximable in general. By contrast, let us show that RNA-PK-Fold(N ), the nearest-neighbor version of the problem, is non-approximable. More precisely, let us show the stronger property that, unless P = N P , there is no polynomial-time algorithm that guarantees to find a matching whose free-energy approximates that of the optimal matching up to a strictly positive factor r(n).

Theorem 4 There exists instances of the nearest-neighbor model such that RNA-PK-Fold(N ) / ∈ AP X.

Proof. Let us briefly outline our proof strategy. We encode any set of numbers X as a string w, having length polynomial on the sum of values in X, and whose matchings are either forbidden (+∞ free-energy), empty (0 free-energy), or have negative energy. Focusing on the latter category, we show that any negative energy matching can be turned, in polynomialtime, into a solution to the 3-PARTITION problem. It follows that any polynomial-time algorithm that guarantees a positive-ratio approximation, thereby producing a matching having negative free-energy anytime such a matching exists, immediately yields a polynomial-time algorithm for the 3-PARTITION problem. The NP-hardness of this problem allows us to conclude on the hardness of approximating RNA-PK-Fold(N * ), within any positive ratio, in the nearest-neighbor energy problem. Let us consider the 3-PARTITION problem, fully defined in Section 3. For any instance X = {x i } 3m i=1 of the problem, let us consider the following RNA sequence:

w = C x 1 AC x 2 A • • • AC x 3m A G K UG K U • • • G K U m times U 2m
Moreover let us consider a nearest-neighbor energy model N * , defined by a function ∆ * N such that:

∆ * N : (A) C C G G -→ -1, ∀i < j, (B) C X Y G -→ -1, ∀i < j, ∀X = C, ∀Y,
(i + 1 and j -1 must both base-pair somewhere, possibly together)

(C) A X Y U -→ -1, ∀i < j, ∀(X, Y ),
(i + 1 and j -1 must both base-pair somewhere, possibly together)

(D) -→ +∞, ∀i < j. i i+1 j-1 j i i+1 j-1 j i i+1 j-1 j

Any other motif

Lemma 5 Let X be a 3-PARTITION instance whose values are bounded by P (n). Then the following statements are equivalent:

-There exists a 3-partition of X into m triplets of equal sum.

-There exists a matching of strictly negative energy over w under N * .

Proof. X is 3-PARTITIONABLE ⇒ ∃m * such that N * w (m * ) < 0: Since X is 3-PARTITIONABLE, then there exists a partition of X into m disjoint triplets ((x a j , x b j , x c j )) m j=1 whose sum are identically K. Consider the matching that pairs each G K block with one of the triplet of blocks C xa j , C

x b j and C xc j , creating nested sequences of base-pairs, and completed with 3 • m (A, U) unconstrained base-pairs over the remaining positions. Clearly, all the positions are involved in a base-pair, and consecutive CC • • • GG are nested as required by energy rule (A). Therefore, any base-pair falls within the scope of energy rules (A), (B) or (C), and the final energy of the matching is

N * w (m * ) = -m • (K + 3) < 0. ∃m * such that N * w (m * ) < 0 ⇒ X is 3-PARTITIONABLE:
Let us start by proving that, within an energy model N * , any matching of w having negative energy is a perfect matching, i.e. every position in the matching is paired. Since N * only allows (C, G) and (A, U) pairs, therefore any valid (finite, negative contribution) base-pair (i, j) must involve a position in the left half of w (C or A) and a position in its right half (G or U), i.e. such that i ≤ m • (K + 3) < j. In order to be valid, (i, j) must also be in a context where i + 1 (resp. j -1) is paired to j ′ ≥ m • (K + 3) (resp. i ′ < m • (K + 3)). The same argument applies to (i + 1, j ′ ) and (i ′ , j -1), and one easily shows by induction that any matching featuring a base-pair (i, j) has infinite energy unless every position in [i, j] is paired. It follows that any matching having negative energy is perfect on some interval Let us finally show that a 3-PARTITION of X can be retrieved from a matching having negative energy. Remind that energy rule (A) forces two consecutive occurrences of C to pair with consecutive occurrences of G. This property extends transitively, and any C x i block in w must therefore be entirely connected to a single G K block. Since a matching of negative energy is perfect, then all the positions in a G K block must be base-paired. Two C x i blocks are not sufficient (x i < K/2) to saturate a G K block, and four blocks would be too large (x i > K/4), violating the constraint that each block must be entirely paired to a single G K block. Therefore a triplet (C xa i , C x b i , C xc i ) blocks is in total interaction with each G K block, and the corresponding values (a i , b i , c i ) constitute a 3-PARTITION of X.

| A = |w [a,b] | U and |w [a,b] | C = |w [a,b] | G . Observe that, since x i < K/2, one has |w [a,b] | A |w [a,b] | C ≤ 1 1+K/2 . Furthermore, if b falls before the final run U 2m , then b < m • (2K + 4) and one has 1 1+K ≤ |w [a,b] | U |w [a,b] | G . It follows that |w [a,b] | A |w [a,b] | C < |w [a,b] | U |w [a,b] | G , i.

⊓ ⊔

From Lemma 5, one knows that the existence of a 3-PARTITION for X can be derived from the existence of a matching of w having negative energy under N * . Now assume there exists a polynomial-time algorithm A that guarantees an r(n) > 0 approximation ratio. Then A would produce a matching m such that N * w (m) = N * w (m * )/r(n), for m * the optimal matching. In particular, A would produce a matching having negative energy anytime such a matching exists. One could then decide, in polynomial time, the 3-partitionability of any set X. Since the decision version of 3-PARTITION is NP-Hard, then there is no such algorithm unless P = N P .

⊓ ⊔

Conclusion/perspectives

We considered the influence of the energy model on the computational complexity of RNA folding with general pseudoknots. In the simplest base-pair model, the problem is exactly equivalent to finding a maximal weighted matching in the graph of compatible positions, and can be solved in Θ(n 3 ) [START_REF] Tabaska | An RNA folding method capable of identifying pseudoknots and base triples[END_REF]. By contrast, it was previously established that the more expressive nearest-neighbor model made the problem NP-Hard [START_REF] Akutsu | Dynamic programming algorithms for RNA secondary structure prediction with pseudoknots[END_REF][START_REF] Lyngsø | RNA pseudoknot prediction in energy-based models[END_REF]. We completed this result by showing that this problem is actually inapproximable within any ratio. Turning to a less expressive -yet realistic -stacking energy model, we have showed that, although NP-hard, the problem could be approximated in polynomial time, at least up to a 1 5 approximation ratio. Quite nicely, a similar approach could be used to refine the computational complexity of RNA-RNA interaction prediction. Already proven NP-complete by Alkan et al [START_REF] Alkan | Rna-rna interaction prediction and antisense rna target search[END_REF], it can be verified that our approximation algorithm achieves the same ratio for RNA-RNA interactions. Furthermore, our NP-hardness and inapproximability results consider bi-partite strings, for which an algorithm for the RNA-RNA interaction problem, suitably parameterized, would yield the same matching as an algorithm for RNA folding with general pseudoknots.

These results show a difference in essence between the nearest-neighbor and the stacking models, which could serve as a starting point for a design of practical (approximation) algorithms for the stacking version of the problem. To that purpose, we plan to complement this study by investigating the existence of a polynomial-time approximation scheme for the problem. Another direction for complementing this study would consider the impact of the energy model on the parameterized-complexity of the problem.

Fig. 2 .

 2 Fig. 2. Typical picture of a standard pseudoknot/matching (Left) and individual contributions of local motifs to the three energy models considered (Right). Dark nodes indicate the supporting base pair for each motif (i.e. (i, j) pairs in our definition below).

Fig. 3 .

 3 Fig.3. Illustrating the reduction: Finding a 3-partition of a set of numbers X (A) is equivalent to finding a matching for wX that produces a maximal number of stacking pairs (C), from which one easily deduces a set of equally summing triplets (D). Such a matching can be represented as a pseudoknotted secondary structure (B).

  [a, b], a ≤ m • (K + 3) < b, and leaves the remaining positions unpaired. Now let us consider which bounds for the interval [a, b] are compatible with a negative energy. Let us denote by w [a,b] the [a, b] factor in a sequence w, and by |w| t the number of occurrences of some letter t in w, then one has |w [a,b]

  e. the matching cannot be perfect on [a, b], and its energy cannot be negative. We are then left to consider the case where m• (2K + 4) ≤ b ≤ |w X |. In such a case, one has |w [a,b] | G = m • Kand one has a = 1. Indeed any greater value a would lead to less thanx i = m • K copies of C, and some G would be left alone. Remark that |w[1,b] | A = 3m, so one must have b = |w|, from which we conclude that any matching having negative energy is perfect, i.e. base-pairs every position.

Acknowledgement

This work was supported by the French ANR MAGNUM ANR 2010 BLAN 0204 grant (YP) and by an INRIA Postdoc Program (SS).