F Lu

M K Iyer

G Parthasarathy

L.-C Wang

K.-T Cheng

K C Chen

An Efficient Sequential SAT Solver With Improved Search Strategies

) a new priority-based search strategy and a flexible sequential search framework which integrates different search strategies, and (3) a decision variable selection heuristic more suitable for solving the sequential problems. We present experimental results to demonstrate that our sequential SAT solver can achieve orders-of-magnitude speedup over Satori.

We plan to release the source code of Seq-SAT along with this paper.

I. Introduction

Boolean SAT finds applications in many areas of circuit design and verification such as Bounded Model Checking [START_REF] Biere | Symbolic Model Checking using SAT Procedures instead of BDDs[END_REF][START_REF] Cimatti | NuSMV 2: An OpenSource Tool for Symbolic Model Checking[END_REF], Unbounded Model Checking [START_REF] Mcmillan | Applying SAT methods in unbounded symbolic model checking[END_REF][START_REF] Kang | SAT-Based Unbounded Symbolic Model Checking[END_REF], Redundancy Identification, Equivalence Checking [START_REF] Tafertshofer | IGRAINE -An Implication-Graph based Engine for Fast Implication, Justification and Propagation[END_REF], Preimage Calculation [START_REF] Li | A novel SAT all-solutions solver for efficient preimage computation[END_REF], etc. State-of-the-art SAT algorithms, as implemented in tools such as ZCHAFF [START_REF] Moskewicz | Chaff: Engineering an efficient SAT solver[END_REF], BERKMIN [START_REF] Goldberg | BerkMin: A fast and robust Sat Solver[END_REF], and C-SAT [START_REF] Lu | A circuit SAT solver with signal correlation guided learning[END_REF] have demonstrated that very hard SAT problems can now be solved in reasonable time. Bounded sequential search using SAT has been shown to be very effective in model checking. However, its major disadvantage is its lack of completeness in general sequential search. A sequential SAT solver Satori was proposed in [START_REF] Iyer | SATORI -A Fast Sequential SAT Engine for Circuits[END_REF], which utilizes combined ATPG and SAT techniques to realize a sequential SAT solver by retaining the efficiency of Boolean SAT and being complete in the search.

Given a sequential circuit, we assume that the circuit follows the Huffman synchronous sequential circuit model as illustrated in Figure 1. Sequential SAT (or sequential justification) is the problem of finding an ordered sequence of input assignments to a sequential circuit, such that a desired objective is satisfied, or proving that no such sequence exists. A desired objective can be a collection of signal value constraints such as that a primary output is 1.

In this paper, we focus on the sequential problems that initial states are given. In sequential SAT, a sequential circuit is conceptually unfolded into multiple copies through backward time-frame expansion (Figure 1). In each time-frame, the circuit becomes combinational and hence, a combinational SAT solver can be applied. In each time-frame, a state element such as a flip-flop is translated into two corresponding signals: a pseudo primary input (PPI) and a pseudo primary output (PPO). The initial state is specified at the PPIs in time-frame 0. The objective is specified at the signals in time-frame n (the last timeframe, where n is unknown before solving the problem). While solving in an intermediate time-frames i (0 < i < n -1), intermediate state solutions are produced at the PPIs and they become intermediate objectives for further justification in the previous time-frames i -1.

Seq-SAT has three major improvements over Satori:

State reduction

As pointed out in the previous work [START_REF] Iyer | SATORI -A Fast Sequential SAT Engine for Circuits[END_REF], minimizing the set of assignments to the state variables, called state reduction, is a critical step for improving the efficiency of sequential search. Seq-SAT performs state reduction not only when an intermediate solution is produced but also when an intermediate state objective is proved unsatisfiable. In the first case, it employs an efficient two-step state reduction algorithm to obtain a smaller state clause. The algorithm follows a similar process as the D-algorithm [START_REF] Abramovici | Digital Systems Testing and Testable Design[END_REF] originally proposed for ATPG and utilizes 3-value simulation. For the second case, a different conflict analysis procedure rather than UIP based conflict analysis [START_REF] Zhang | Efficient conflict driven learning in a Boolean satisfiability solver[END_REF] is applied to derive a smaller state clause.

Flexible sequential search framework Sequential search can

be carried out in two strategies: depth-first or breadth-first. In the depth-first strategy, a SAT algorithm expands the search to a new timeframe whenever a solution state is identified for the current timeframe. If the backward justification fails in timeframe i -1 and no state solution can be found, it then backtracks to timeframe i and select another state solution of timeframe i to continue the justification. In the breadth-first strategy, the state solutions in the current timeframe i are exhausted before expanding into a new timeframe i -1.

In this paper, we provide a flexible sequential search framework in which these two search strategies can be integrated. We introduce a priority-based search strategy that can combine both depthfirst search and breadth-first search and demonstrate its efficiency over other search strategies.

Decision variable selection heuristic

In SAT, the idea that make decisions to satisfy recently deduced clauses [START_REF] Moskewicz | Chaff: Engineering an efficient SAT solver[END_REF][START_REF] Goldberg | BerkMin: A fast and robust Sat Solver[END_REF] has been proven to be very effective for solving hard combinational problems. However, when applying a combinational solver for sequential SAT, the run time of combinational SAT solving is not necessarily the most important factor to be optimized. In this paper, we propose a decision variable selection heuristic which is more suitable for solving sequential problems.

The rest of the paper is organized as the following. Section II illustrates the use of combinational SAT in sequential SAT. In Section III, we discuss the importance of state reduction and describe our approaches to this problem. Section IV presents our sequential search framework. In Section V, we present our priority-based search strategy. Section VI discusses the impact of decision ordering in combinational SAT on the overall sequential SAT efficiency. Section VII summarizes experimental results to demonstrate the superiority of our current sequential SAT solver. Section VIII concludes the paper.

II. Combinational and sequential SAT

Given a finite set of variables, V , over the set of Boolean values B ∈ {0, 1}, a literal, l/l is an instance of a variable or its complement,

v/¬v ∈ V . A clause c i , is a disjunction of literals (l 1 ∨ l 2 ... ∨ l n ). A formula f , is a conjunction of clauses c 1 ∧ c 2 ... ∧ c n .
A clause is considered as a set of literals, and a formula as a set of clauses.

An assignment A satisfies a formula f if f (A) = 1. An assignment is called maximal when every variable in V receives a value assignment. Given a sequential circuit, a state clause is a clause consisting only of state variables.

In a Boolean Satisfiability (SAT) problem (or here call it combinational SAT to differentiate from sequential SAT), a formula f is given and the problem is to find an assignment A to satisfy f or prove that no such an assignment exists. SAT has attracted tremendous research effort in recent years, resulting in the developments of various efficient SAT solver packages such as [START_REF] Moskewicz | Chaff: Engineering an efficient SAT solver[END_REF][START_REF] Goldberg | BerkMin: A fast and robust Sat Solver[END_REF][START_REF] Lu | A circuit SAT solver with signal correlation guided learning[END_REF][START_REF] Zhang | SATO: An Efficient Propositional Prover[END_REF][START_REF] Marques-Silva | GRASP: A Search Algorithm for Propositional Satisfiability[END_REF]. Through backward timeframe expansion (Figure 1), a sequential SAT problem can be translated into a sequence of combinational SAT problems. To illustrate the usage of state clauses in sequential SAT, Figure 2 depicts a simple example circuit with three primary inputs a, b, c, one primary output f , and three state elements x, y, z. The initial state is (x = 1, y = 0, z = 1). Suppose the SAT objective is to satisfy f = 1.

II-A. Sequential SAT and state clauses

Starting from time-frame n where n is unknown, the circuit is translated into a combinational copy with state variables duplicated as PPOs and PPIs. This is illustrated as (1) in the Figure . Since this represents a combinational SAT problem, a Boolean SAT solver can be applied.

Suppose after the combinational SAT solving, a solution (a = 1, b = 1, c = 0, PPI x = 0, PPI y = 1, PPI z = 0) for satisfying f = 1 (step (2))is identified. Notice that all inputs have values assigned in this solution, even though it might not be necessary to have all of them assigned in order to satisfy f = 1. This phenomenon is due to the use of a combinational SAT solver as mentioned earlier.

The assignments at PPIs implies a solution state (x = 0, y = 1, z = 0) that doesn't contain the initial state.

Instead of justifying the state (x = 0, y = 1, z = 0) in time-frame n-1, we can examine it first to check whether we can change as many of the assignments to don't-cares while still satisfying the objective f = 1. Suppose after analysis, we determine that z = 0 is unnecessary, it can be removed from the solution state (x = 0, y = 1, z = 0) and a solution state (x = 0, y = 1) is derived. We call this step solution state reduction since it tries to remove the unnecessary state variable assignments from a solution. The solution state (x = 0, y = 1) doesn't contain the initial state either, so it still needs to be justified in timeframe n -1. Backward time-frame expansion can be achieved by adding a state objective (PPO x = 0, PPO y = 1) to the combinational copy of the circuit. Also, a state clause (x + y ) is generated to prevent reaching state solutions contained by the state sub-space (x = 0, y = 1) in time-frame n -1. The new combinational SAT instance is then passed to the Boolean SAT solver.

Suppose in time-frame n -1, no solution can be found for state objective (PPO x = 0, PPO y = 1). Then, we backtrack to time-frame n to find another solution. In a way, we have proved that from state (x = 0, y = 1), there exists no solution. Therefore, there is no need to remove the state clause (x + y ). However, at this point it is beneficial to perform further analysis to determine if both PPO x = 0 and PPO y = 1 are involved in the conflict that indicates state objective (PPO x = 0, PPO y = 1) is unsatisfiable. Suppose after conflict analysis, It is discovered that only PPO x = 0 is involved in the conflict. A state clause "(x)" is then added. The future combinational SAT instances will then have the state clauses "(x + y )(x)" included, that record the non-solution state sub-spaces previously identified. This is illustrated in (4) of the Figure . We call this step objective state reduction since it tries to derive a smaller state clause when a state objective becomes unsatisfiable.

The solving continues until either one of the following two conditions is reached:

1. After solution state reduction, a solution is found whose state part contains the initial state. For example, a solution with the state part (x = 1, z = 1) contains the initial state (x = 1, y = 0, z = 1). In this case, a solution for the sequential SAT problem is found. Note that a solution without any assignments to PPIs contains any initial state.

2. If backtracked to time-frame n, the initial objective f = 1 cannot be satisfied under the constraints imposed by all the state clauses added, then the original problem is unsatisfiable. This is equivalent to say: if any objective (including the initial objective f = 1 and all intermediate objectives) cannot be satisfied under the constraints imposed by all the state clauses added, then the original problem is unsatisfiable.

The above example illustrates several important concepts in the design of our current sequential SAT solver.

• The state reduction involves finding smaller state clauses in order to more effectively prune the search space. There are two types of state reduction. One is solution state reduction, the other is objective state reduction. • The use of state clauses serves two purposes: (1) to record those state sub-spaces that have been explored, and (2) to record those state sub-spaces containing no solution. The first is to prevent the search from entering a state justification loop. • Although conceptually the search follows the backward timeframe expansion, the above example demonstrates that in the implementation, explicit time-frame expansion is not necessary.

In other words, a sequential SAT solver needs only one copy of the circuit. Moreover, there is no need to memorize the number of time-frames being expanded. Later we will discuss our implementation of the solver and demonstrate how this can be achieved by using an objective list. In this framework, the keys are to decide how to select the next objective from the list for solving and to decide how many solutions to be found for an objective each time. The details will be given in Section IV.

In the following, we will describe the detail of Seq-SAT. To facilitate the description, we define the following two terms: (1) a frame objective is an objective to be satisfied which is passed to the combinational SAT solver. A frame objective can be either the initial objective or a state objective. (2) a frame solution is an assignment at the PIs and PPIs, which satisfies a given frame objective.

III. State reduction for state space pruning

There are two types of state reduction: solution state reduction (SSR) and objective state reduction (OSR). SSR is performed when a frame solution is identified by the combinational SAT solver. OSR is applied when a state objective becomes unsatisfiable. The details are discussed in the following.

III-A. Solution state reduction

When a frame solution is identified by the combinational SAT solver, all PIs and PPIs have assigned values. Our goal is to find and un-assign the unnecessary assignments at the PPIs. With fewer assigned PPIs, a smaller state clause can be created, which prunes a larger state sub-space.

Given an initial frame solution, to derive one containing a minimal number of state variables with value assignment is an NP-complete problem. Therefore, we only derive a locally minimized 3-value solution. A frame solution s of a frame objective o is a locally minimized 3-value solution if and only if un-assign any state variable of s (changing its value from a 0/1 to X), the resulting assignments cannot satisfy o by 3-value simulation.

One straight forward heuristic to derive a locally minimized 3value solution is to employ 3-value simulation. That is, for all assigned state variables of a solution, un-assign one of them at a time to check if its assignment is necessary to satisfy the frame objective. If not, the assignment is removed from the solution, this procedure continues until each assigned state variable has been checked once, then the remaining assignments form a 3-value minimized solution.

A frame solution reported by a SAT solver usually contains many unnecessary assignments at state variables. Directly apply 3-value simulation could be computationally expensive. For example, removing N unnecessary state assignments needs at least N times of simulation even if parallel simulation is used. To reduce the number of simulation runs and improve efficiency, we adopt a two-step state reduction approach. In this approach, a modified D-algorithm [START_REF] Abramovici | Digital Systems Testing and Testable Design[END_REF] is first applied to obtain a minimized solution. Given a frame objective O and its frame solution S, the modified D-algorithm justify O based on S, that is, when the D-algorithm makes a decision at a signal, the value assigned to the signal has to be consistent with its value in the solution S. Hence, the D-algorithm is used as a trace procedure, not as a search procedure. Based on the minimized solution from the first step, 3-value simulation is then employed to derive a 3-value minimized solution. This approach is based on the assumption that the D-algorithm can usually find a solution containing much fewer assignments than that of a SAT solver. Experiments show our assumption is valid in general. The details of the approach are omitted due to page limit. Table I shows experimental results to demonstrate the impact of solution state reduction. All our experiments were run on P4 2GHz Linux machines with 1 GB memory.

In these experiments, a sequence of sequential SAT problem instances are created for each circuit. For each circuit and each of its primary outputs f , the first SAT objective is to satisfy f = 1 and then we would try to satisfy f = 0. This is a typical primary output toggle experiment. The results shown for each circuit are for solving all toggle SAT instances for the circuit. We assume the initial state is that all state elements have a 0 value. Note, this kind of problems may contain very difficult cases even for ISCAS 89 circuits.

The two "#" columns show the numbers after state reduction with modified D-algorithm (step 1) and then, with three value simulation (step 2), respectively. The two "% reduced" columns show the corresponding percentage of reduction. The two run-time columns show the run times based on using step 1 and then, using step 2. We note that if no state reduction is performed, then the sequential SAT might not finish solving all instances within a reasonable time for each of these examples.

From the results, we can see that the modified D-algorithm could reduce most of the unnecessary assignments to the state variables, so the number of 3-value simulation runs can be greatly reduced. Although step 2 does not help much in terms of the total percentage of assignment reduction, it does help further improve the overall performance.

III-B. Objective state reduction

As a state so is added to the objective list, a state clause sc is added to prevent the state from being reached again. When the state objective so becomes unsatisfiable, instead of using the same state clause sc, a smaller state clause might be added to prune a larger search space. This can be achieved by additional conflict analysis when a state objective is proved unsatisfiable. We explain the idea by an example, suppose the state objective so is (PPO r = 1, PPO s = 1, PPO t = 0, PPO x = 0, PPO y = 1), the implication graph when sc is proved unsatisfiable (which means there is a conflict that can not be resolved by backtrack) is shown as in Figure 3. From Figure 3, if tracing from the conflict points to the PPO state variables, we can deduce that (PPO r = 1, PPO s = 1, PPO t = 0) is unsatisfiable, so the state clause (r + s + t), which is smaller than (r + s + t + x + y ), can be added to the circuit. However, we need to be careful that the state clauses added in this way shouldn't exclude the search space that contains the initial state. Otherwise a satisfiable problem erroneously becomes unsatisfiable since the initial state couldn't be reached anymore after adding such a state clause. Using the above example, if the values of r, s and t in the initial state are 1, 1 and 0 respectively, state clause (r + s + t) would exclude the initial state. To avoid this, we can just select a state variable from the state objective whose objective value is not its initial value and add its corresponding literal to the state clause. For example, if the value of x in the initial state is 1 in the above example, we can add x to the state clause (r + s + t) to form state clause (r + s + t + x ) which can be added to the circuit. This can always be done since state objectives don't contain the initial state. That is, each state objective must have at least one variable whose objective value is not same as its value in the initial state, otherwise, the problem should has been proved satisfiable.

This reduction is valid even if state clauses that were used for preventing justification loop are involved in the conflict graph.

We also note that this reduction is not related to identify the subset of a CNF formula that is sufficient for unsatisfiability, instead, it just tries to find a conflict clause that contains only literals of the variables in a state objective.

Table II gives the experimental results comparing the performance of our sequential solver with and without OSR. These are the same toggle experiments as those in Table I. From the experiment results, it can be seen that the performance is greatly improved for half of the testcases by OSR, and only for one testcase, it causes a small overhead.

IV. Flexible sequential search framework

In our sequential SAT design, state clauses are accumulated through the solving process. The sequential solving process consists of a sequence of combinational solving tasks based on a given frame objective. At the beginning, the frame objective is the initial objective. As the solving proceeds, many frame solutions become frame objectives. These frame objectives are stored in the objective list.

The use of state clauses and the objective list offer much flexibility in the choice of the search strategy. Different search strategies can be implemented by using different criteria in selecting the next frame objective for solving, and in deciding how many solutions to be obtained at a time for a frame objective. For example, for the depth-first search, the objective list is used as a stack where the next selected objective is the one most recently produced. For the breadth-first search, the objective list is used as a queue. In addition, in the breadth-first search, it demands the solver to find all solutions for the selected objective. After that, the objective is removed from the objective list. This framework further enables us to explore various search strategies between these two extreme ones which will be described later in this section.

The overall framework of our sequential SAT solver is described in Algorithm IV.1. In the description, we assume that each frame objective produces at most one frame solution at a time. This assumption is only for purpose of easier explanation. The algorithm can be easily modified to serve the case that multiple solutions are produced at a time. The procedure PPO state conflict analysis() does objective state reduction, and state reduction() does solution state reduction. The procedure select a frame objective() selects a frame objective to be solved next, upon different implementation, various search strategies can be realized, and the details will be given in the next section.

A frame objective can be removed from the objective list only if it is proven unsatisfiable by the combinational SAT. If it is satisfiable, the frame objective stays in the objective list. If the objective list becomes empty during the sequential search, it means that the solver has exhausted all the state objectives and, therefore, proved that the initial objective is unsatisfiable.

Note that during each step of combinational SAT, conflict clauses are also accumulated through the combinational SAT solving process. When the sequential solving switches from one frame objective to another, these conflict clauses stay. Hence, in the sequential solving process, the conflict clauses generated by the combinational SAT are also accumulated. Our experience indicates that although these conflict clauses does help speed up the combinational SAT solving, for sequential SAT, the state clauses are the dominating source for the efficiency improvement of the overall sequential search.

Algorithm IV.1: SEQUENTIALSOLVER(C,ob j,s 0 ) comment: C is the circuit with PPIs and PPOs expanded comment: ob j is the initial objective comment: s 0 is the initial state comment: FO is the objective list

FO ← {ob j}; while (FO = / 0) do                                      f ob j ← select a frame objective(FO); f sol ← combinational solve a frame objective(C, f ob j ); if ( f sol = NULL)    clause ← PPO state conflict analysis(C, f ob j ); add state clause(C,clause); FO ← FO -{ f ob j }; else              state cube ← state reduction(C, f ob j , f sol ); if (s 0 ∈ state cube) return (SAT ); else    clause ← convert

V. On the scheduling of frame objectives

In this section, we introduce a new search strategy, called the priority-based search strategy, which can be easily implemented in the search framework described above. We further compare it with the depth-first and breadth-first search strategies and demonstrate its superiority.

Given a state objective, we define its size as the total number of state variables minus the number of state variables with an binary objective value (i.e. not an unknown X) in the state objective. We use size(ob j) to denote the size of ob j. Note that a state objective with a larger size corresponds a larger state sub-space. Hence, if there are multiple state objectives in the objective list, it seems intuitive to choose the largest state objective as the next objective for solving as it corresponds to exploring the solution with the largest state sub-space, and its solutions, if any, tend to contain fewer state assignments, which means smaller state clauses are likely to be generated to prune a larger search space. We further define the count of an objective ob j as the number of times that ob j has been selected for solving by the combinational solver. If we use count to guide the selection of objectives and choose the one with the largest count as the next objective, it will be more like the breadth-first search strategy. On the other hands, if we choose the one with the smallest count, It will be more like the depth-first search. Now we define priority of a frame objective ob j, pri(ob j), as follows:

pri(ob j) = size(ob j)count(ob j)/10

(1) In our implementation of the priority-based search strategy, we select the objective with the largest value of pri. If there is a tie under this selection criterion, then the most recently generated one will be selected. For each selected objective, at most one solution is obtained at a time. Through extensive experiments, we observe that combining size and count in this way maximizing the performance, in general, for both satisfiable and unsatisfiable cases.

With the objective list, depth-first and breadth-first search strategies can be implemented easily. In depth-first search, the list behaves like a stack. Hence, the selected objective is the one most recently produced. In breadth-first search, the list behaves like a queue. However, in breadth-first search, solving an objective means to finish producing all possible solutions. Then, the objective is removed from the objective list.

Table III gives the experimental results comparing different search strategies. These are the same toggle experiments as those in Table I. The "priority" column is the same as the last "times" column in Table I. As it can be seen, the new search strategy significantly outperforms the depth-first and the breadth-first strategies.

VI. On the decision variable selection heuristics in the combinational SAT solver

For pure combinational SAT, both VISDS [START_REF] Moskewicz | Chaff: Engineering an efficient SAT solver[END_REF] decision variable selection and clause based decision variable selection [START_REF] Goldberg | BerkMin: A fast and robust Sat Solver[END_REF] have been proven effective. However, when applying a combinational SAT solver for sequential SAT, the runtime of combinational SAT solving is not the only critical factor. In each run of combinational SAT solving during the sequential SAT, the smaller the necessary assignments are made at the PPIs, the more efficiently the state space can be pruned. Therefore, the decision variable selection heuristics in combinational SAT also need modifications. To demonstrate this point, we conducted experiments to compare different decision variable selection heuristics in combinational SAT.

Decision heuristic 1: This is the same as the VISDS decision variable selection [START_REF] Moskewicz | Chaff: Engineering an efficient SAT solver[END_REF].

Decision heuristic 2: First, the decision variable is selected from the most recently generated conflict clause that is not yet satisfied [START_REF] Goldberg | BerkMin: A fast and robust Sat Solver[END_REF]. If all the generated conflict clauses are satisfied, the decision variable is selected based on the J-nodes [START_REF] Ganai | Combining strengths of circuit-based and CNF-based algorithms for a high-performance SAT solver[END_REF][START_REF] Lu | A circuit SAT solver with signal correlation guided learning[END_REF] with the highest topological order (closest to the POs and the PPOs). This decision heuristic is used in the original combinational SAT solver C-SAT [START_REF] Lu | A circuit SAT solver with signal correlation guided learning[END_REF]. Decision heuristic 3: First, the decision variable is selected from the most recently generated conflict clause that is not yet satisfied as in Heuristic 2 above. If all the generated conflict clauses are satisfied, the decision variable is selected from PIs and PPIs where PIs are selected before PPIs This is the default heuristic used in our final implementation.

The reasons we adopt Heuristic 3 in Seq-SAT are summarized as follows:

(1). Using input variables as decision points results in better performance for easy combinational problems in general. For a sequential problem that can be transformed into a series of relatively easy combinational problems, this heuristic works very well. (2). For hard problems, as the conflict clauses tend to accumulate fast during search, the decisions will gradually be dominated by the conflict clauses. Therefore, the performance will not degrade much. (3). Selecting PIs prior to PPIs as decision variables tends to result in solutions with fewer state variables with value assignments [START_REF] Kang | SAT-Based Unbounded Symbolic Model Checking[END_REF]. Table IV shows the results of different heuristics. The experiments for "s-" circuits were the same PO toggle experiments circuit as explained before. The experiments for "c-" circuits were the combinational equivalence checking problem using the miter circuit of two identical ISCAS 85 circuits.

It can be seen that the decision heuristic 3 works well for sequential circuits, even for large combinational circuits, its performance is not bad.

VII. Additional experimental results

Table V compares our solver with NuSMV [START_REF] Cimatti | NuSMV 2: An OpenSource Tool for Symbolic Model Checking[END_REF] bounded model checking function for satisfiable cases. We have shown in the previous experiments how each idea works, Table VI compares the overall performance of Seq-SAT with Satori.

In Table V, the "Property" column gives the LTL formula for each testcase, and "length" columns give the length of the witness vectors obtained by Seq-SAT and NuSMV respectively. For NuSMV, the version we used is NuSMV 2.1.2-zchaff in which the SAT solver zchaff is used for bounded model checking. The initial states used in these examples are all 0. From Table V, it can be seen that, even though the witness vectors derived by Seq-SAT are in general longer than those by NuSMV, Seq-SAT clearly outperforms NuSMV in terms of CPU runtime.

In Satori, the depth-first search and the breadth-first search are implemented separately. So in Table VI, we compare Seq-SAT with Satori under these two search strategies for the same toggle experiments as those in Table I It can be seen that Seq-SAT achieves significant performance improvement for almost all cases. It is also interesting to note that the results of the depth-first search in Satori are quite different from those of depth-first search implemented in Seq-SAT (see Section V Table III). This indicates that the depth-first search is quite unstable (especially for satisfiable cases) and easy falls into bad search area.

Figure 1 :

 1 Figure 1: Backward time-frame expansion

Figure 2 :

 2 Figure 2: Sequential SAT with state clauses

Figure 3 :

 3 Figure 3: An implication graph when a state objective is proved unsatisfi able

  There will be three possible outcomes reported by the solver: (1) o is unsatisfiable, (2) at least one state in the obtained solutions contains the initial state s 0 , and (3) no state in the obtained solutions contains the initial state s 0 . For case (1), we can remove o from FO. For case (2), the search stops and the sequential problem is proven satisfiable. For case (3), the state part of each solution is added to FO as a new state objective and a state clause is added to C to prevent the state be explored and found again in any future solutions. If finally, the objective list FO becomes empty, which means all objectives have been proven to be unsatisfiable, then the search stops and the sequential problem is proven unsatisfiable.

• The above example demonstrates the use of depth-first search strategy. However, by using state clauses, more flexible search strategies are feasible if the intermediate state objectives are recorded in a list. The idea of the sequential search framework can be described as follows: Given a sequential problem with an objective obj and an initial state s 0 , C is the one time-frame combinational copy of the sequential circuit where each state variable is expanded to a PPI variable and a PPO variable. Suppose we use an objective list FO to store the intermediate objectives, and initially it only contains obj. Each time, an objective o is selected from FO and given to the combinational circuit solver.

The circuit solver then solves o on C.

TABLE I :

 I RESULTS OF SOLUTION STATE REDUCTION

	Circuit Org #	Step 1	resulting	Step 2	resulting
			#	%reduced times(sec)	#	%reduced times(sec)
	s526	91134 70981 22.1%	1	70545	0.6%	1
	s1423 237344 73986 69.1%	3	73298	0.4%	1
	s5378 181487 89499 50.6%	9	87450	2.3%	2
	s35932 3382	3142	7.1%	2	3142	0%	3
	s38417 2008744 1054999 47.4%	280	1052562 0.2%	100
	s38584 869942 375770 56.8%	128	371707 1.1%	110

TABLE II :

 II EFFECTS OF OBJECTIVE STATE REDUCTION

		Without OSR	With OSR
	Circuit runtimes(sec)	stats	runtimes(sec)	stats
	s526	2	12/0/0	1	12/0/0
	s1423	1	10/0/0	1	10/0/0
	s5378	402	89/5/4	3	90/8/0
	s13207	1924	165/58/19	1345	166/63/13
	s15850	2002	100/54/20	522	100/69/5
	s35932	3	640/0/0	3	640/0/0
	s38417	286	210/0/2	100	212/0/0
	s38584	106	498/58/0	111	498/58/0
	**abort time for each problem is 100 sec.
	stats: x/y/z are #SAT/#UNSAT/#Abort

TABLE III :

 III CPU RUNTIMES (IN SECONDS) ON SEARCH STRATEGIES

		depth-fi rst	breadth-fi rst	priority
	Circuit runtimes	stats	runtimes	stats	runtimes	stats
		(sec)		(sec)		(sec)
	s526	9	12/0/0	2	12/0/0	1	12/0/0
	s1423	301	7/0/3	103	9/0/1	1	10/0/0
	s5378	1504	75/8/15	511	89/5/4	3	90/8/0
	s13207 2159 167/54/21 1367 166/63/13 1345 166/63/13
	s15850 1401 100/60/14 1569 100/59/15	522	100/69/5
	s35932 32020 320/0/320	3	640/0/0	3	640/0/0
	s38417 2324	189/0/23	2031	192/0/20	100	212/0/0
	s38584 8515 414/57/85 2006 482/57/17	111	498/58/0
	**abort time for each problem is 100 sec.
	stats: x/y/z are #SAT/#UNSAT/#Abort

TABLE IV :

 IV CPU RUNTIMES ON DIFFERENT DECISION HEURISTICS

		heuristic 1	heuristic 2	heuristic 3
	Circuit runtimes	stats	runtimes	stats	runtimes	stats
		(sec)		(sec)		(sec)
	s526	1	12/0/0	3	12/0/0	1	12/0/0
	s1423	205	8/0/2	1	10/0/0	1	10/0/0
	s5378	313	90/7/1	109	90/7/1	3	90/8/0
	s13207	1436	166/62/14	1374	166/63/13	1345	166/63/13
	s15850	511	100/69/5	514	100/69/5	522	100/69/5
	s35932	3	640/0/0	3	640/0/0	3	640/0/0
	s38417	194	211/0/1	253	210/0/2	100	212/0/0
	s38584	41	498/58/0	300	496/58/2	111	498/58/0
	c3540	40	0/1/0	13	0/1/0	20	0/1/0
	c5315	21	0/1/0	16	0/1/0	8	0/1/0
	**abort time for each problem is 100 sec.
	stats: x/y/z are #SAT/#UNSAT/#Abort

TABLE V :

 V COMPARISON TO NUSMV BOUNDED MODEL CHECKING

			Seq-SAT		NuSMV
	Circuit	Property	runtimes(sec) length runtimes(sec) length
	s526	G(g214 = 0)	1	133	15	33
	s13207 G(g594 = 0)	1	5	2	5
	s13207 G(g785 = 0)	7	1 1	34	10
	s38417 G(g5549 = 0)	9	815	> 3600	-
	s38417 G(g16399 = 0)	14	7	55	7
	s38584 G(g11678 = 0)	6	3 0	831	12
	s38584 G(g29212 = 0)	48	144	> 3600	-
	**abort time for each problem is 3600s.

  . Table VII summarizes the number of the abort cases shown in Table VI. Column "#Satori abort" shows the number of cases aborted in both Satori depth-first search and Satori breadth-first search. Column "#SSANA" gives the number of cases aborted in Seq-SAT but not aborted in either Satori depth-first search or Satori breadth-first search.

TABLE VI :

 VI COMPARISON TO SATORI

		Satori depth-fi rst	Satori breadth-fi rst	Seq-SAT
	Circuit runtimes	stats	runtimes	stats	runtimes	stats
		(sec)		(sec)		(sec)
	s526	405	8/0/4	5	12/0/0	1	12/0/0
	s1423	400	6/0/4	106	9/0/1	1	10/0/0
	s5378	4100	53/4/41	1455	80/4/14	3	90/8/0
	s13207	4010	154/48/40	1836	166/59/17	1345	166/63/13
	s15850	2302	99/52/23	2214	100/52/22	522	100/69/5
	s35932	3200	608/0/32	1	640/0/0	3	640/0/0
	s38417	805	204/0/8	2257	191/0/21	100	212/0/0
	s38584	4810	451/57/48	2228	477/58/21	111	498/58/0
	**abort time for each problem is 100 sec.
	stats: x/y/z are #SAT/#UNSAT/#Abort
		TABLE VII: ABORTED CASES ANALYSIS
		Circuit #Satori abort #Seq-SAT abort #SSANA
		s526	0		0	0
		s1423	1		0	0
		s5378	14		0	0
		s13207	17		13	0
		s15850	22		5	0
		s35932	0		0	0
		s38417	5		0	0
		s38584	11		0	0

Proceedings of the Design, Automation and Test in Europe Conference and Exhibition (DATE'05) 1530-1591/05 $ 20.00 IEEE

Proceedings of the Design, Automation and Test in Europe Conference and Exhibition (DATE'05)

VIII. Conclusion

In this paper, we present an efficient sequential SAT solver Seq-SAT. The implementation of Seq-SAT employs four new ideas: (1) better state reduction algorithms, (2) a more flexible search framework for accommodating and integrating different search strategies, (3) the priority-based search strategy, and (4) a modified decision variable selection heuristic for the underlying combinational circuit SAT solver. With these new ideas and efficient data structures and implementation, Seq-SAT achieves very significant speedup over Satori.

We will release the source code of Seq-SAT along with the publication of this paper.


Sonny Tarbouriech

Wael Suleiman

On Bisection Continuous Collision Checking Method: Spherical Joints and Minimum Distance to Obstacles

 to efficiently handle the case of spherical and two revolute joints, this kind of joints is very common in modern robotic systems. The new formulations provide more tight motion bounds, thus increase the success rate of checking collision-free paths. We also propose an extension to get the minimum distance to obstacles along a path, this information is primordial as it allows samplingbased motion planning techniques to sort collision-free paths according to their minimum clearance.

We have integrated our implementation into a samplingbased motion planning technique and validated it through simulation and on the real Baxter research robot. The experiments revealed that the method not only does not miss any collision between the robot and the obstacles, but also the minimum distance extension provides the path with the maximum clearance at no additional computational cost.

I. INTRODUCTION

Collision checking is an essential step in motion planning as it ensures the path to be collision-free, and it has been a very active research topic during the last two decades. However, this is also one of the most computationally expensive operations in the process of finding a clear path, especially in a complex environment.

The main challenge relies on determining whether the continuous path between two states in the Configurationspace (C-space) is in collision or not. In sampling-based algorithms, a commonly used method consists in testing individual configurations by linear interpolation between the states. The principal cumbersome of such an approach is that a compromise must be made when selecting the number of samples. Decreasing it leads to a loss of accuracy and could yield to a collision omission. On the other hand, increasing it could as well increase the computational time for collision checking, which is practically inefficient strategy in situations where the robot is far from obstacles.

In the context of unstructured environments where the safety of the robot, and possibly of humans interacting with it, is a priority, a computationally efficient method for exact collision checking is primordial.

The guarantee of a collision-free path has been introduced by Quinlan in [START_REF] Quinlan | Real-time modification of collision-free paths[END_REF] and is obtained by using overlapping free bubbles in the C-space along a segment. This concept was then integrated in several motion planning strategies [START_REF] Vahrenkamp | Efficient motion planning for humanoid robots using lazy collision checking and enlarged robot models[END_REF], [START_REF] Vahrenkamp | Rdt+: A parameter-free algorithm for exact motion planning[END_REF]. {Sonny.Tarbouriech@USherbrooke.ca, Wael.Suleiman@USherbrooke.ca} The concept of Continuous Collision detection (CCD) has also received a lot of attention as it guarantees a path to be collision-free. A CCD algorithm for articulated models has been proposed in [START_REF] Redon | Fast continuous collision detection for articulated models[END_REF], it is able to compute the time of contact and contact locations. A fast CCD algorithm for articulated models using Taylor models and temporal culling has been proposed in [START_REF] Zhang | Continuous collision detection for articulated models using taylor models and temporal culling[END_REF], and another fast and exact CCD algorithm which uses properties of the Bernstein basis and Bézier curves has been presented in [START_REF] Tang | Fast and exact continuous collision detection with bernstein sign classification[END_REF]. Moreover, an algorithm that performs CCD for polygon-soup models undergoing rigid motion at highly interactive rates is presented in [START_REF] Tang | C2A: Controlled conservative advancement for continuous collision detection of polygonal models[END_REF]. CCD algorithms have been also used in local planning in samplingbased motion planning [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF], [START_REF] Redon | Practical local planning in the contact space[END_REF], [START_REF] Schwarzer | Exact collision checking of robot paths[END_REF].

In this paper, we are interested in the category of CCD for motion planning, mainly the CCD method proposed in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF] and [START_REF] Schwarzer | Exact collision checking of robot paths[END_REF], this method is overviewed in Section II.

Our contributions are twofold:

• Adapting the CCD method in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF] to the case of two revolute and spherical joints; hence proving more tight motion bounds. This contribution is presented in Sections III-A and III-B. • Proposing an extension to [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF] to get information about the minimum distance to obstacles along a collision-free path. The detail of this contribution is given in Section IV.

We have also validated our algorithm on a Baxter research robot through simulation and real experiments. The experimental results are analyzed in Section V.

II. ADAPTIVE BISECTION CONTINUOUS COLLISION CHECKING METHOD

The main idea behind this method is to establish a sufficient condition of collision-free by computing the geometric path of rigid bodies in the workspace. Lemma 1 in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF] states that a sufficient condition to guarantee that two rigid objects, A 1 and A 2 , do not collide at any configuration q located on the path segment π, which is joining two configurations q a and q b , is to verify the following inequality:

λ 1 (q a , q b ) + λ 2 (q a , q b ) < η 12 (q a ) + η 12 (q b ) (1)
where η 12 (q i ) is the minimum distance between objects A 1 and A 2 for a given configuration q i , and λ i (q a , q b ) refers to the maximum Euclidean displacement of all the points in an object A i along the path segment π.

In the case of a robotic arm, we will focus on the two following situations (represented in Fig. 1):

• Self-collision detection: Both objects A 1 and A 2 are distinct links of the robot. • Detection of collisions with the environment: A 1 is a link of the robot and A 2 is an obstacle in the world.

As it is an off-line analysis, we assume that objects in the environment are static. Consequently, in this case, λ 2 (q a , q b ) = 0 in (1).

To make a complete collision check for a given path segment, we have to evaluate all pairs of elements which could potentially collide using [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF]. If the test succeeds for a pair, we ensure that those two objects do not collide and we move on to the next pair. If the inequality is not verified, we cannot affirm that the path is collision-free. However, that does not mean a collision necessarily exists, but more tests are needed. To this end, we add an intermediate state on the path segment π. This new configuration divides the initial segment in two parts. If a collision is detected for this state, the motion is then invalid. Otherwise, we evaluate in the same manner each subsegment independently for the same two objects. If both subsegments satisfy the inequality, the path segment is collision-free for the given pair of elements, if not, we repeat the process for the subsegments that failed the test. More details about the method are given in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF]. The procedure of applying this method to a multi-DOF robotic arm is also presented in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF].

III. DISPLACEMENTS IN C-SPACE AND WORKSPACE

To apply the method in Section II, the term λ i (q a , q b ) in (1) should be computed to verify the inequality constraint. Computing the exact displacement of a multi-DOF robotic arm between two configurations could be a complex process, and therefore a time-consuming operation. However, an upper bound of this distance can be quickly figured out by only considering the displacement range of each joint and the dependency between the links movements.

The original method in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF] considers the case of a single link L of length L, which is controlled by a revolute joint. In this case, an elementary joint variation dθ leads to a traveling distance for the link of

dλ = Ldθ (2)
Fig. 2 illustrates the total displacement of the link between configurations q a and q b . The point P ∈ L which travels the longest distance is located at the end of the link and moves from P A to P B . The total length s is:

λ = θ b θ a Ldθ = L(θ b -θ a ) (3)

A. Two Revolute Joints Case

The case where two revolute joints control the same link, as illustrated in Fig. 3, is particularly interesting for modern redundant robotic arms, and it is the case for the Baxter research robot. 3. Displacement of a link controlled by two revolute joints By using the original formula proposed in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF], the upper limit of the displacement would be:

θ 1 L θ 2 P A P B Fig.
λ * = L |θ b 1 -θ a 1 | + |θ b 2 -θ a 2 | (4)
It is clear that this over-estimation could increase the probability that a segment fails the collision-free test in (1), and as a result a sub-segmentation is required.

To compute a more tighter upper-bound, let us start by computing an elementary displacement of a point P situated at the end of the link:

dλ = L (dθ 2 ) 2 + sin 2 θ 2 (dθ 1 ) 2 (5)
The total traveled distance can be defined by:

λ = θ b 1 θ a 1 θ b 2 θ a 2 L (dθ 2 ) 2 + sin 2 θ 2 (dθ 1 ) 2 (6)
By linearly interpolating between the configurations q a and q b , we obtain dθ 2 = βdθ 1 , where β is a constant. Let us suppose that sin 2 θ * 2 = max θ a 2 ≤θ2≤θ b 2 sin 2 θ 2 , thus the integral ( 6) can be upper bounded by a simple integral:

λ ≤ λ = θ b 1 θ a 1 L β 2 + sin 2 θ * 2 dθ 1 = L β 2 + sin 2 θ * 2 θ b 1 -θ a 1 = L θ b 2 -θ a 2 2 + sin 2 θ * 2 θ b 1 -θ a 1 2 (7)
It is clear that the new upper limit in ( 7) is more tighter than the one in ( 4), an exemple with the Baxter research robot to compare those upper limits is given in Appendix. For this linkage, the upper bounds λi (i = 1, 2) can be written as follows:

B. Spherical Joints Case

λ1 =L 1 θ b 2 -θ a 2 2 + sin 2 θ * 2 θ b 1 -θ a 1 2 λ2 =(L 1 + L 2 ) θ b 2 -θ a 2 2 + sin 2 θ * 2 θ b 1 -θ a 1 2 + L 2 (θ b 3 -θ a 3 ) + L 2 θ b 5 -θ a 5 2 + sin 2 θ * 5 θ b 4 -θ a 4 2 (8
) where sin θ * × , θ a × and θ b × are defined in the same way as in Section III-A.

Note that λ 1 is the same displacement defined in [START_REF] Tang | Fast and exact continuous collision detection with bernstein sign classification[END_REF], this is because θ 3 does not affect the workspace displacement of the first link. On the other hand, θ 3 affects the displacement of the second link. This result can be easily extended to the case of an articulated model with N spherical joints.

IV. MINIMUM OBSTACLE DISTANCE COMPUTATION

In [START_REF] Schwarzer | Exact collision checking of robot paths[END_REF], it has been proven that a minimum clearance ζ between objects A 1 and A 2 at any configuration q between q a and q b is guaranteed if:

λ 1 (q a , q b ) + λ 2 (q a , q b ) < η 12 (q a ) + η 12 (q b ) -2ζ (9)
This feature is interesting, however it suffers from several disadvantages, for instance:

1) As the minimum clearance is fixed, the planner might ignore a narrow passage, which would be the only feasible solution.

2) The criterion in [START_REF] Redon | Practical local planning in the contact space[END_REF] does not allow to compare two clear paths. However, in motion planning we are usually interested in sorting admissible paths according to their clearance, this is particularly interesting for safety issue, either for the robot or for a human interacting with the robot. For a path where two objects A 1 and A 2 move between two configurations q a and q b , we define the estimated clearance as:

δ = η 12 (q a ) + η 12 (q b ) -λ 1 (q a , q b ) -λ 2 (q a , q b ) 2 (10)
If A 2 is a fixed obstacle, λ 2 (q a , q b ) = 0 and (10) becomes:

δ = η 12 (q a ) + η 12 (q b ) -λ 1 (q a , q b ) 2 = dist(q a , q b , A 1 , A 2 ) (11)
It is worth pointing out that the value of δ can be either negative or positive:

• δ ≤ 0 implies that the motion cannot be guaranteed to be collision-free, this is because (1) is not satisfied. A sub-segmentation is necessary. • δ > 0, then δ is the minimum distance between A 1 and A 2 along the path.

λ 1 q a q b η 12 η 12 2|δ| (a)
The potential maximum traveled distance by the link is higher than the sum of minimum distance from the object at each configuration (δ ≤ 0)

q a q b η 12 η 12 λ 1 2|δ| (b)
The two circles fully cover the link displacement and a margin distance from the obstacle is guaranteed (δ > 0) Fig. 5. Illustration of the relations between distances in the workspace.

A schematic example that illustrates the above two cases is given in Fig. 5. For more details about this representation, refer to [START_REF] Schwarzer | Exact collision checking of robot paths[END_REF].

The procedure to compute the minimum clearance along a path segment and sorting collision-free segment paths according to their clearance is given in Algorithm 1. Note that each element of the structure segment refers to a specific Algorithm 1 Minimum obstacle distance along a path procedure CLEARANCE(q a , q b ) priority queue Q segments sorted by increasing δ for i ← 1 : nb links do for j ← 1 : nb obstacles do δ ← dist(q a , q b , A i , A j ) Equation ( 11) if (η ij (q a ) ≤ 0)||(η ij (q b ) ≤ 0) then return COLLISION else Q.Add Segment(q a , q b , i, j, δ) end if end for end for repeat segment s ← Q.T op() Q.P op T op() q new ← Interpolate(s.q a , s.q b ) δ 1 = dist(s.q a , q new , s.link, s.obst) δ 2 = dist(q new , s.q b , s.link, s.obst) if η s.link,s.obst (q new ) ≤ 0 then return COLLISION else Q.Add Segment(s.q a , s.q new , s.link, s.obst, δ 1 ) Q.Add Segment(s.q new , s.q b , s.link, s.obst, δ 2 ) end if until Q.T op().λ ≤ λ is defined in (1) if Q.T op().δ ≤ 0 then return COLLISION else return Q.T op().δ end if end procedure pair of link/obstacle evaluated between two states and is used to store the corresponding distance information. Parameter can be defined as the maximum admissible error in the distance estimation. It is a positive user-defined constant that affects the performances of the algorithm: decreasing it improves the returned distance estimation accuracy whereas increasing it reduces the required computational burden to generate the estimation, the impact of on the accuracy and the computational complexity is shown in Section V. Table I.

As it has been mentioned in Section II, if the inequality (1) is not satisfied, a sub-segmentation by adding an intermediate state is necessary. However, the location of that state can have a significant impact on the algorithm convergence rate. Many covering strategies have been proposed in [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF], however we selected the strategy that samples toward the part of the path that is most likely to be in collision. This strategy yields to a minimum number of segments evaluation, and can be described as follows :

• When δ < 0, the new sampled state corresponds to the middle of the uncovered length (Fig. 6(a)) • When δ ≥ 0, the new sampled state is the projection of the intersection of the two circles onto the path segment (Fig. 6

(b))

We refer to this strategy as the Interpolate(q a , q b ) procedure in Algorithm 1. Lemma. The estimated and exact distances to obstacles satisfy the following inequality:

δ exa - 2 ≤ δ ≤ δ exa ( 12
)
where δ exa and δ are, respectively, the exact and estimated minimum distances between two objects A 1 and A 2 , where A 1 moves from configuration q a to q b and A 2 is a fixed obstacle.

Proof: Let q worst be the configuration, between q a and q b , where the distance is minimal.

From [START_REF] Sucan | Moveit![END_REF], we can define the error between the exact value and our estimation as:

∆ = δ exa -δ = η 12 (q worst ) - η 12 (q a ) + η 12 (q b ) -λ 1 (q a , q b ) 2 (13)
The proof of Equation ( 12) can now be constructed as follows:

• ∆ = 0 is the minimum error: By definition, [START_REF] Sucan | Moveit![END_REF] gives an underestimation of the minimum clearance along a path segment [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF]. That means δ exa ≥ δ est is always verified and thus ∆ ≥ 0. The value ∆ = 0 can be obtained in at least one case: when q a = q b . That means that, in the best case, our method gives no error on the estimation.

• ∆ = 2 is the maximum error:

To get the maximum value of ∆, (13) can be rewritten as follows:

∆ = η 12 (q worst ) - η 12 (q a ) + η 12 (q b ) 2 + λ 1 (q a , q b ) 2 (14)
We now have two independent parts in the expression of ∆. As η 12 and λ 1 are always positive, the maximum value of ∆ is obtained when the two following conditions are fulfilled:

η 12 (q a ) + η 12 (q b ) is minimal λ 1 (q a , q b ) is maximal As the minimum clearance is obtained at configuration q worst , that means η 12 (q * ) ≥ η 12 (q worst ) where * is a or b. As defined in Algorithm 1, the stopping criterion of the procedure is λ ≤ . That means the maximum value that λ 1 (q a , q b ) can take is . When replacing these terms in ( 14), we get ∆ = 2 . This is the maximum error on the estimated clearance that is obtained in the worst case.

Hence, the condition δ exa -2 ≤ δ ≤ δ exa is always satisfied.

V. EXPERIMENTAL RESULTS

The efficiency and performance of the described method have been tested in simulation and on the real robot Baxter from Rethink Robotics Inc. Baxter is equipped with two 7-DOF manipulator arms. In this study, we only focus on the control of one arm but our method can be easily extended to also consider both arms simultaneously. Our algorithm has been implemented using the platform MoveIt! [START_REF] Sucan | Moveit![END_REF] that integrates the motion planning library OMPL [START_REF] S ¸ucan | The Open Motion Planning Library[END_REF] and communicates via ROS (Robotic Operating System) API interface [START_REF] Quigley | ROS: an open-source robot operating system[END_REF].

A. Distance Approximation Function

Even though the continuous collision checking method tends to limit the number of distance computations, those operations remain widely called during the planning process. Moreover, computing the exact minimum distance between two complex and concave objects is a computationally expensive operation. Moveit! uses the Flexible Collision Library (FCL) [START_REF] Pan | FCL: A general purpose library for collision and proximity queries[END_REF] as its primary collision checking library. Besides collisions check, FCL can also compute the distance between two non-overlapping objects.

To accelerate the execution of our algorithm, we developed a computationally efficient distance computation function that gives a reasonable approximation of the distance between a link of the robot and an obstacle. This function is based on simplifying the environment representation by transforming each obstacle into its corresponding oriented bounding box (OBB). The Baxter robot is also approximated by its collision geometry model from the Unified Robot Description Format (URDF) file. Thus, the arms are modeled as cylinders and boxes. Each part of the arm is then subdivided into an optimal number of spheres that encompass the initial volume. Finally, the minimum distance between an obstacle and a link becomes the minimum distance between the OBB and the spheres and, hence, can be easily computed.

B. Numerical Experiments

The first validation phase of our approach was made through simulated scenarios using the Gazebo simulator [START_REF] Koenig | Design and use paradigms for gazebo, an open-source multi-robot simulator[END_REF] as an interface for the Baxter robot. All experiments are performed on Intel(R) Core(TM) i5-3470 CPU @ 3.20 GHz PC with 8 GB RAM.

Two scenarios have been considered in this study:

• Scenario 1: The right arm is queried to reach successively two particular configurations for the joints. These states are located on the sides of two balloons tethered to the ground by thin rods as shown in Fig. 7(a). • Scenario 2: The right arm's end-effector is planned to reach specific poses (position and orientation) passing through a wire grid leaving only narrow free spaces see Fig. 7(b). The conventional method of uniformly discretizing a path to check if it is collision-free is called hereafter the discrete method and our implementation as continuous method.

We implemented the continuous and discrete collision checking methods into the classical RRT* and their performances in finding the minimum distance from obstacles along a path is compared in in Table I. The tests have been carried out in the two scenarios to point out the environmental influence. Each single test gathers 100 motion planning tasks, and the allowed time for task planning is set to 20 seconds for both implementations. The standard distance function (FCL) has been used in both cases, and for different values of . For the discrete case, is the sampling rate in C-space. It refers to the length of the longest segment that does not require state validation and is defined as a fraction of the space maximum extent. For the continuous method, is the maximum acceptable error in the distance evaluation between two objects.

It is worth pointing out that, for a single arm of Baxter robot, the maximum extent in the C-space is approximately 10π (the sum of the joint ranges). Thus, for a standard value = 5%, a segment will be considered as collision free if the total corresponding joint variation is under π/2. In the worst case, this variation can produce a displacement of more than 1, 7m. Conversely, if we want to ensure a maximum error of 10mm in the workspace, as done with the continuous method, we have to set a value of = 0.03%. This is not an efficient solution as the computational time will tremendously increase.

The distance obtained by the discrete method, with an extremely low value of (0.01%), is taken as a reference. We indicate by ∆ the error between the reference and the estimated distance.

The sign of ∆ plays a major role with respect to the safety issue. As it is always positive for the continuous method, we ensure that no collision could occur.

The computation time is only indicative because the accuracy of both methods is not comparable. Note that the use of the distance approximation function can tremendously improve the computational time compared to FCL (with an average gain of 74.6 times in Sc.1 and 5.5 times in Sc.2). However, this strategy underestimates the distance between the objects and might add a non-negligible error on the evaluation (7.1mm for Sc.1 and 52.5mm for Sc.2 on average). This is a good option when planning in relatively uncluttered areas because it is very fast to compute. On the other hand, that approximation might neglect narrow passages as they could be considered in collision. Table II presents the failure rates as a result of an undetected collision. We chose a fixed value = 10mm for the continuous collision checker to enable the robot to operate in a complex environment, like the one in Scenario 2.

As anticipated, the continuous method did not omit any collision and enabled the planner to achieve all tasks. On the other hand, the discrete method tends to find inadmissible solutions due to undetected collisions.

We have also validated the experiment in Scenario 1 on the real Baxter robot as shown in Fig. 8. The whole experiment can be seen in the accompanying video.

VI. CONCLUSION

This paper presented a new formulation to bound the motions of spherical and two revolute joints. This formulation is integrated into the bisection collision checking method [START_REF] Schwarzer | Adaptive dynamic collision checking for single and multiple articulated robots in complex environments[END_REF] in order to provide tighter motion bounds than the original formulation. An extension to get the minimum distance to obstacles is also presented, this extension can be directly integrated into any sampling-based motion planning algorithm to sort collision-free paths according to their clearance.

Our experiments show that the computational time for the new implementation of bisection collision checking method is comparable or faster than the conventional discrete method, most importantly it never miss a potential collision with the environment. For those reasons, it can be an appropriate candidate to remplace the conventional discrete method in environments with thin objects.

Future work will focus on integrating our implementation of the bisection collision checking method into a motion planning algorithm to:

• Exploit the feature of minimum distance to obstacles to maximize human safety in the context of human-robot cooperation.

• Deal with bi-manipulation tasks by integrating the second arm of the Baxter research robot.

APPENDIX

To compare the upper limits in [START_REF] Tang | Fast and exact continuous collision detection with bernstein sign classification[END_REF] and (4), we consider the shoulder's two revolute joints (joints 1 and 2 in Fig. 9) of the Baxter robot arm, and the following initial and final configurations:

θ a 1 = - π 2 , θ b 1 = π 3 θ a 2 = - π 3 , θ b 2 = π 3 (15)
The arm length starting from the axe of joint 1 is L = 0.9679 m.

In this cas, the upper limits are: where λ and λ * are the upper limit computed by ( 7) and ( 4) respectively.

The above results point out the significant difference between the two upper limits on the workspace displacement, and that our new formulation provides a tighter bound on the displacement than the original one.

  Fig. 1. Example: Two types of collision analysis for a 3-DOF robotic arm

Fig. 2 .

 2 Fig. 2. Displacement of a link controlled by one joint

1 ✓ 2 ✓ 3 L 2 Fig. 4 .

 12324 Fig. 4. Displacement of two links controlled by two spherical joints

when δ ≥ 0 Fig. 6 .

 06 Fig. 6. Illustration of the covering strategy

Fig. 7 .

 7 Fig. 7. Motions that connect near states in the C-space could create a large displacement in the workspace. Discrete collision checker does not detect collisions (represented by red spheres) with the thin rod.

λ = 3

 3 .4997 m λ * = 5.3213 m ∆λ λ * -λ = 1.8216 m

Fig. 8 . 2 Joint 1 LFig. 9 .

 8219 Fig.8. Snapshots of the Baxter robot executing a motion planning in an unstable environment using our implementation of the continuous collision checking method.

ACKNOWLEDGMENT

This research is supported by the Natural Sciences and Engineering Research Council of Canada (NSERC) under the grant RGPIN-2012-419406.


Evgenia Kartsaki

Bruno Cessac

Gerrit Hilgen

Evelyne Sernagor

How specific classes of retinal cells contribute to vision: A Computational Model

Vision begins with the photoreceptors converting light from the visual scene into electrical signals, compressing our visual world into a code of action potentials sent to the brain by the retinal ganglion cells (RGCs). A human retina contains almost 1 million RGCs and each of these cells interprets different features of the visual scene (shape, motion, color, etc.). It is all these parallel streams of information received by the brain, that eventually lead to visual perception. Currently, there exist over 30 RGCs subtypes based on: Ø common anatomical features, Ø functional properties, Ø common gene expression.

Contemporary questions:

✓What role does each RGC subtype play in vision? ✓How is vision impaired if one of these subtypes is inactivated?

☛ We propose a novel approach combining for the first time pharmacogenetics, electrophysiology, morphology, behavior and mathematical modelling in order to selectively inactivate specific RGCs types and decipher their role in vision, both at the single cell and population level.

dV dt =( -g L (V-V L ) -g Na m 3 h (V-V Na ) -g K n 4 (V-V K ) -g CNO (V-V K ) 𝑑𝑚 𝑑𝑡 = 𝑎

IntroductionPharmacogenetics:

  Modify neuronal activity noninvasively and reversibly by using the Designer Receptors Exclusively Activated by Designer Drugs (DREADDs) technique. DREADD-expressing RGCs' activity can be reduced or silenced when the synthetic designer drug clozapine N-oxide (CNO) is added to the solution. Light stimulation and electrophysiological recordings: Record light responses from the RGC layer at pan retinal level in vitro with the APS MEA in control and CNO conditions Morphology: Visualize cells and their dendritic trees by imaging under confocal microscopy post-recording. Behavioral testing: Perform behavioral tests in control and CNO conditions ☛ How does CNO impact the response of RGCs to visual stimuli? RGCs response to light stimuli is characterized by their spatio-temporal receptive field (RF). Spatial profile: the region of the visual field in which light stimuli evoke responses in a RGC. The simplest RF is organized into a center-surround structure, responding oppositely to light. Temporal profile: the temporal course of this response. RF estimation using PRANAS Estimated receptive field of an OFF RGC (A) In control condition (B) In CNO condition Estimated receptive field of an ON RGC (A) In control condition (B) In CNO condition Spatial profile Temporal profile Develop a generalized RGC model Model single cell dynamics with Hodgkin-Huxley conductances using the NEURON simulation platform in order to reproduce what we observe in the experiments. Disclaimer: All the experimental work is done at Pr. Sernagor's lab Active Pixel Sensor MEA Chip featuring 4,096 electrodes(42 μm spacing) arranged in a 64x64 configuration, covering an active area of 7.12 mm 2 . From Maccione et al., J Physiol 2014 Pan-retinal recordings in control (left) and CNO (right) conditions Methods ☛ Use mathematical modelling The difference of Gaussians function can model the spatial RF ☛ Use numerical simulations PRANAS: A platform for retinal analysis and simulation, developed by the Biovision team at Inria. MACULAR: A platform for large scale simulations of the retina in pathological conditions (currently under co-development by the Biovision team and the CED engineering team at Inria). NEURON: A simulation platform for modelling individual and networks of neurons, developed at Yale and Duke.

C

  9 𝑉 1 -𝑚 -𝑚 𝑉 𝑚 𝑑ℎ 𝑑𝑡 = 𝑎 = 𝑉 1 -ℎ -𝛽 = 𝑉 ℎ 𝑑𝑛 𝑑𝑡 = 𝑎 @ 𝑉 1 -𝑛 -𝛽 @ 𝑉 𝑛 Source: https://benchling.com/pub/dreadds From Prusky et al., Iovs 2004 The center-surround mechanism is temporally biphasic Conclusions Ø We reproduced the strong effect of CNO on the temporal part, as in the experiments. Next steps: -Relate the CNO conductance to the CNO concentration -Fit the experimental data, to tune the model's parameters Ø There isn't a clear effect on the spatial part yet. -There might be, if amacrine cells are also affected by CNO Ø Mimic experimental setup with MACULAR. is reduced RGC's activity is silenced To study the effect of CNO on the behavior of the cell, we added a leak potassium channel depended on [CNO], represented by the extra current g CNO .

This research is funded by the Leverhulme Trust in the context of the project "A novel approach to functional classification of retinal ganglion cells", involving the Institute of Neuroscience, Newcastle University and the Biovision team, INRIA SAM.