Accelerating Parametric Probabilistic Verification

Accelerating Parametric Probabilistic Verification ⋆ Nils Jansen1 , Florian Corzilius1 , Matthias Volk1 , Ralf Wimmer2 , Erika Ábrahám1 , Joost-Pieter Katoen1 , and Bernd Becker2 1 RWTH Aachen University, Germany {nils.jansen | corzilius | volk | abraham | katoen}@cs.rwth-aachen.de 2 Albert-Ludwigs-University Freiburg, Germany {wimmer | becker}@informatik.uni-freiburg.de Abstract. We present a novel method for computing reachability probabilities of parametric discrete-time Markov chains whose transition probabilities are fractions of polynomials over a set of parameters. Our algorithm is based on two key ingredients: a graph decomposition into strongly connected subgraphs combined with a novel factorization strategy for polynomials. Experimental evaluations show that combining these approaches leads to a speed-up of up to several orders of magnitude in comparison to existing approaches. 1 Introduction Discrete-time Markov chains (DTMCs) are a widely used modeling formalism for systems exhibiting probabilistic behavior. Their applicability ranges from distributed computing to security and systems biology. Efficient algorithms exist to compute measures like: “What is the probability that our communication protocol terminates successfully if messages are lost with probability 0.05?”. However, often actual system parameters like costs, faultiness, reliability and so on are not given explicitly. For the design of systems incorporating random behavior, this might even not be possible at an early design stage. In modelbased performance analysis, the research field of fitting [1], where—intuitively— probability distributions are generated from experimental measurements, mirrors the difficulties in obtaining such concrete values. This calls for treating probabilities as parameters and motivates to consider parametric DTMCs, PDTMCs for short, where transition probabilities are (rational) functions in terms of the system’s parameters. Using these functions one can, e. g., find appropriate values of the parameters such that certain properties are satisfied or analyze the sensitivity of reachability probabilities to small changes in the parameters. Computing reachability probabilities for standard DTMCs is ⋆ This work was partly supported by the German Research Council (DFG) as part of the research project CEBug (AB 461/1-1), the Research Training Group AlgoSyn (1298), the EU FP7-project MoVeS, the FP7-IRSES project MEALS and by the Excellence Initiative of the German federal and state government. typically done by solving a linear equation system using iterative methods. This is not feasible for PDTMCs. Instead, approaches based on state elimination have been proposed [2,3]. The idea is to replace a state and its incident transitions with direct transitions from its predecessor to its successor states. In this way, one eliminates all states except for the initial and target states of a system. The result is a rational function describing the probability of reaching a set of target states, depending on the values of the parameters. The efficiency of such elimination methods strongly depends on the order in which states are eliminated and on the representation of the rational functions. Related work. The idea of constructing a regular expression representing the DTMC originates from Daws [2]. He uses state elimination to generate regular expressions describing the paths to the target states of the system. Hahn et al. [3] apply this idea to PDTMCs to obtain rational functions for reachability and expected reward properties. They improve the efficiency of the construction by heuristics for the transformation of finite automata to regular expressions [4] to guide the elimination of states. Additionally, they reduce the polynomials to simplify the rational functions. These ideas have been extended to Markov decision processes [5]. The main problem there is that the reachability probabilities depend on the chosen resolution of the nondeterminism. When maximizing or minimizing these probabilities, the optimal resolution generally depends on the values of the parameters. Their algorithms are implemented in PARAM [6], the—to the best of our knowledge—only available tool for computing reachability probabilities of PDTMCs. This paper can be seen as a continuation of [2,3]. Several authors have considered the related problem of parameter synthesis: for which parameter instances does a given (LTL or PCTL) formula hold? To mention a few, Han et al. [7] considered this problem for timed reachability in continuous-time Markov chains, Pugelli et al. [8] for Markov decision processes (MDPs), and Benedikt et al. [9] for ω-regular properties of interval Markov chains. Contributions of the paper. In this paper we improve the computation of reachability probabilities for PDTMCs [2,3] in two important ways. We consider a state elimination strategy based on a recursive graph decomposition of the PDTMC into strongly connected subgraphs and give a novel method to efficiently factorize polynomials. Although presented in the context of parametric Markov chains, this constitutes a generic method for representing and manipulating rational functions and is also suited for other applications as well. The experiments show that the combination of both techniques yields a speed-up of more than one order of magnitude compared to [3]. 2 Preliminaries Definition 1 (Discrete-time Markov chain). A discrete-time Markov chain (DTMC) is a tuple D = (S, I, P ) with a non-empty finite set S of states, an P initial distribution I : S → [0, 1] ⊆ R with Ps∈S I(s) = 1, and a transition probability matrix P : S × S → [0, 1] ⊆ R with s′ ∈S P (s, s′ ) = 1 for all s ∈ S. The states SI = {sI ∈ S | I(sI ) > 0} are called initial states. A transition leads from a state s ∈ S to a state s′ ∈ S iff P (s, s′ ) > 0. The set of successor states of s ∈ S is succ(s) = {s′ ∈ S | P (s, s′ ) > 0}. A path of D is a finite sequence π = s0 s1 . . . sn of states si ∈ S such that P (si , si+1 ) > 0 for all i < n. Paths D is the set of all paths of D. We denote the set of paths starting at s ∈ S by Paths D (s). Paths D (s, t) denotes the set of paths starting at s and ending at t. We generalize this to sets S ′ , S ′′ ⊆ S of states by Paths D (S ′ , S ′′ ) = S S D ′ ′′ D s′′ ∈S ′′ Paths (s , s ). A state t is reachable from s iff Paths (s, t) 6= ∅. s′ ∈S ′ D D The probability measure Pr fin (π) for paths π ∈ Paths is given by Pr D fin (π =
Qn−1 s0 . . . sn ) = I(s0 ) · i=0 P (si , si+1 ). Note that for two paths π1 , π2 ∈ Paths D D D it holds that Pr D fin ({π1 , π2 }) = Pr fin (π1 ) + Pr fin (π2 ) if no path is a prefix of the other one, i. e., if they are independent. For a set R ⊆ P stochastically D (R) = (π) with R′ = {π ∈ R | ∀π ′ ∈ Paths D we define Pr D Pr fin fin π∈R′ ′ R. π is not a proper prefix of π}. For more details we refer, e. g., to [10]. For a DTMC D = (S, I, P ) and a subset of states K ⊆ S we define the set of input states of K by Inp(K) = {s ∈ K | I(s) > 0 ∨ ∃s′ ∈ S \ K. P (s′ , s) > 0}, i. e., the states inside K that have an incoming transition from outside K. Analogously, we define the set of output states of K by Out(K) = {s ∈ S \ K | ∃s′ ∈ K. P (s′ , s) > 0}, i. e., the states outside K that have an incoming transition from a state inside K. The set of inner states of K is given by K \ Inp(K). We call a state set S ′ ⊆ S absorbing iff there is a state s′ ∈ S ′ from which no state outside S ′ is reachable in D, i. e., iff Paths D ({s′ }, S r S ′ ) = ∅. A state s ∈ S is absorbing if {s} is absorbing. A set S ′ ⊆ S induces a strongly connected subgraph (SCS) of D iff for all s, t ∈ S ′ there is a path from s to t visiting only states from S ′ . A strongly connected component (SCC) of D is a maximal (w. r. t. ⊆) SCS of S. If Out(S ′ ) = ∅ holds for an SCC S ′ , S ′ is called a bottom SCC. The probability of eventually reaching a bottom SCC in a finite DTMC is always 1 [10, Chap. 10.1]. We only consider probabilistic reachability properties, i. e., the probability
to D eventually reach a set T ⊆ S of target states, formally: Pr D (S , T ) . It Paths I fin is well-known that this suffices for checking arbitrary ω-regular properties, see [10, Chap. 10.3] for the details. 2.1 Parametric Markov Chains To add parameters to DTMCs, we follow [6] by allowing arbitrary rational functions defining probability distributions. Definition 2 (Polynomial and rational function). Let V = {x1 , . . . , xn } be a finite set of variables with domain R. A polynomial g over V is a sum of monomials, which are products of variables in V and a coefficient in Z: e e g = a1 · x11,1 · . . . · xen1,n + · · · + am · x1m,1 · . . . · xenm,n , where ei,j ∈ N0 = N ∪ {0} and ai ∈ Z for all 1 ≤ i ≤ m and 1 ≤ j ≤ n. Z[x1 , . . . , xn ] denotes the set of polynomials over V = {x1 , . . . , xn }. A rational 0.2 0.8 s2 0.4 s1 s3 q 0.2 0.2 1 1−q s4 S1.1 0.4 0.5 S1.2 S1.2.1 s9 1 s7 p 0.3 S1 1 0.2 0.8 s6 s5 s8 1−p Fig. 1. Example PDTMC and its SCC decomposition function over V is a quotient f = gg12 of two polynomials g1 , g2 over V with  g2 = 6 03 . We use FV = gg21 | g1 , g2 ∈ Z[x1 , . . . , xn ] ∧ g2 6= 0 to denote the set of rational functions over V . Definition 3 (PDTMC). A parametric discrete-time Markov chain (PDTMC) is a tuple M = (S, V , I, P ) with a finite set of states S, a finite set of parameters V = {x1 , . . . , xn } with domain R, an initial distribution I : S → FV , and a parametric transition probability matrix P : S × S → FV . As we are applying graph-based algorithms, we need the underlying graph of a  (P)DTMC M, which is given by GM = (S, DP ) where DP = (s, s′ ) ∈ S × S | P (s, s′ ) 6= 0 . Using an evaluation, all or some of the parameters occurring in the rational functions of a PDTMC can be instantiated. Definition 4 (Evaluated PDTMC). An evaluation u of V is a function u : V → R. The evaluation g[u] of a polynomial g ∈ Z[x1 , . . . , xn ] under u : V → R substitutes each x ∈ V by u(x), using the standard semantics for + and ·. For [u] ∈ R if g2 [u] 6= 0. f = gg21 ∈ FV we define f [u] = gg21 [u] For a PDTMC M = (S, V , I, P ), the evaluated PDTMC is the DTMC D = (Su , Iu , Pu ) given by Su = S and for all s, s′ ∈ Su , Pu (s, s′ ) = P (s, s′ )[u] and Iu (s) = I(s)[u] if the evaluations are defined and 0 otherwise. An evaluation u substitutes each parameter by a real number. This induces a probability measure on the evaluated PDTMC under the following conditions. Definition 5 (Well-defined evaluation). An evaluation u is well-defined for PDTMC M = (S, V , I, P ) if for the evaluated PDTMC D = (Su , Iu , Pu ) it holds that P – Pu : Su × Su → [0, 1] with ∀s ∈ Su : s′ ∈Su Pu (s, s′ ) = 1, and 3 g2 6= 0 means that g2 cannot be simplified to 0. – Iu : Su → [0, 1] with P s∈Su Iu (s) = 1. A well-defined evaluation u is called graph preserving, if it holds that ∀s, s′ ∈ S : P (s, s′ ) 6= 0 =⇒ P (s, s′ )[u] > 0. Note that P (s, s′ )[u] > 0 implies that no division by 0 will occur. This will be ensured during the model checking algorithm. Evaluation u is required to be graph preserving, i. e., GM = GMu . This is necessary as by altering the graph, states might become unreachable which can change reachability probabilities. Definition 6. Given a PDTMC M = (S, V , I, P ) and a set of absorbing target states T ⊆ S, the parametric probabilistic model checking problem is to find for each initial state sI ∈ SI and each t ∈ T a rational function fsI ,t ∈ FV such that for all graph-preserving evaluations u : V → R and the evaluated PDTMC D = (Su , Iu , Pu ) it holds that fsI ,t [u] = Pr Mu (Paths(sI , t)). 3 Parametric Model Checking by SCC Decomposition In this section we present our algorithmic approach to apply model checking to PDTMCs. In the following, we assume every PDTMC M = (S, V , I, P ) to have only bottom SCCs consisting of one state, i. e., absorbing states, which will be the target states. For each initial state sI ∈ SI and each target state t ∈ T we compute a rational function fsI ,t over the set of parameters V which describes the probability of reaching t from sI as in [3]. A similar method was introduced in [11] for the non-parametric case. 3.1 PDTMC Abstraction The basic concept of our model checking approach is to replace a non-absorbing subset of states K ⊆ S and its transitions inside a PDTMC M by transitions directly leading from the input states Inp(K) of K to the output states Out(K) of K. These transitions have the probabilities of all paths visiting only states of K. This concept is illustrated in Figure 2: In Figure 2(a), an arbitrary, non-absorbing set of states K has one input state sI and two output states s1out , s2out . The abstraction in Figure 2(c) hides every state of K except for sI ; all transitions are directly leading to the output states. As we need a probability measure for arbitrary subsets of states, we first define sub-PDTMCs induced by such subsets. Definition 7 (Induced PDTMC). Given a PDTMC M = (S, V , I, P ) and a non-absorbing subset K ⊆ S of states, the induced PDTMC over K and M is given by MK = (S K , I K , V K , P K ) with S K = K ∪ Out(K), V K = V , ∀s ∈ S K . I K (s)6=0 ⇐⇒ s ∈ Inp(K), and   P (s, s′ ), if s ∈ K, s′ ∈ S K K ′ 1, if s = s′ ∈ Out(K) P (s, s ) =  0, otherwise . s2out K sI K s1out sI s1out (a) Initial PDTMC K s1out sI s2out (b) Abstraction of K with abstract loop s2out (c) Abstraction of K Fig. 2. Concept of PDTMC abstraction Intuitively, all incoming and outgoing transitions are preserved for inner states of K while the output states are made absorbing. We allow an arbitrary input distribution I with the only constraint that I(s)6=0 iff s is an input state of K. Example 1. Consider the PDTMC M in Figure 1 and the state set K = {s7 , s8 }. The induced PDTMC MK = (S K , I K , V K , P K ) over K and M shown in Figure 3(a) has output states Out(K) = {s5 , s6 , s9 } and input states Inp(K) = {s7 }. For our abstraction we take into account all finite paths that do not intermediately return to the initial state. In Figure 2(b), there are abstract transitions leading to the output states together with a self-loop on the initial state. The outgoing transitions describe all paths that do not visit the input state again, while the self-loop describes all paths that return to the input state. These paths build the set of all paths that add to the probability of finally reaching one of the output states. Note that inside a non-absorbing set of states, the probability of reaching the set of all output states is 1. Figure 2(c) shows the final abstraction where the probability of the self-loop is taken into account in determining the transition probabilities of the outgoing transitions. Formally, we define the probability of all finite paths that start in a state s and finally reach a state s′ without returning to s beforehand. This includes paths that both start and end in s. Definition 8. Assume a PDTMC M = (S, V , I, P ), a non-absorbing state s ∈ S and a state s′ ∈ S. The path abstraction of s and s′ is given by
M M ′ pM (s, s′ ) | si 6= s ∧ si 6= s′ , 0 < i < n} . abs (s, s ) = Pr fin {π = s0 . . . sn ∈ Paths Using this we are now ready to define the abstraction of a PDTMC M with respect to initial states and target states. The probabilities are the total reachability probabilities between these states. Let us first consider an example. Example 2. Consider the PDTMC M′ = (S ′ , I ′ , P ′ , V ′ ) of Figure 3(a) and let the set of target states T ′ = {s5 , s6 , s9 } correspond to the absorbing states of M′ . The s6 s6 1 fs7 ,s6 1 0.5 0.2 s7 s5 s7 1 fs7 ,s7 fs7 ,s5 s8 1−p fs7 ,s9 s9 (a) Induced PDTMC 1 fˆs7 ,s6 1 p 0.3 s6 1 s5 s7 1 fˆs7 ,s5 s5 1 fˆs7 ,s9 1 s9 s9 (b) Abstracted PDTMC (c) Scaled functions Fig. 3. PDTMC Abstraction ′ ′ ′ ′ ′ abstract PDTMC M′abs = (Sabs ) has states Sabs = {s5 , s6 , s7 , s9 } , Iabs , Pabs , Vabs and edges from s7 to all other states. The first abstraction step according to ′ the path abstraction pM abs as in Definition 8 is depicted in Figure 3(b). The rational functions describing the probabilities of all finite paths that either leave K without visiting state s7 again or starting and ending in s7 are: ′ fs7 ,s6 = pM abs (s7 , s6 ) = 0.5 ′ fs7 ,s9 = pM abs (s7 , s9 ) = 0.3 · (1 − p) fs7 ,s5 = pM abs (s7 , s5 ) = 0.2 fs7 ,s7 = pM abs (s7 , s7 ) = 0.3 · p ′ ′ The total probability of reaching the output states is given by paths which first use the loop on s7 arbitrarily many times (including zero times) and then take a transition to an output state. For example, using the geometric series, the probability of the set of paths leading from s7 to s5 is given by ∞ X (fs7 ,s7 )i · fs7 ,s5 = i=0 1 · fs7 ,s5 1 − fs7 ,s7 As the probability of finally reaching the set of absorbing states in M′ is 1, we can directly scale the probabilities of the outgoing edges such that their sum is equal to 1. This is achieved by dividing each outgoing probability by the sum of all outgoing probabilities, fout = 0.2 + 0.5 + 0.3 · (1 − p). The abstract and scaled PDTMC is depicted in Figure 3(c) with the probabilities given by fˆs7 ,s5 = 0.2 /fout fˆs7 ,s9 = (0.3 · (1 − p)) /fout fˆs7 ,s6 = 0.5 /fout We now define the final abstraction formally. Definition 9 (Abstract PDTMC). For a PDTMC M = (S, V , I, P ) and a set of absorbing states T ⊆ S, the abstract PDTMC Mabs = (Sabs , Vabs , Iabs , Pabs ) is given by Sabs = {s ∈ S | I(s) 6= 0 ∨ s ∈ T }, we define Iabs (s) = I(s) and  ′ pM  abs (s, s )  P   M ′′ s′′ ∈T pabs (s, s ) Pabs (s, s′ ) = 1    0 Vabs = V , and for all s, s′ ∈ Sabs if I(s) > 0 ∧ s′ ∈ T if s = s′ ∈ T otherwise . Theorem 1. For a PDTMC M = (S, V , I, P ) and its abstraction Mabs = (Sabs , Iabs , Vabs , Pabs ) according to Definition 9 it holds for all initial states sI ∈ SI and all absorbing states t ∈ T that M Mabs abs Pr M (sI , t)) = Pr M (sI , t)). fin (Paths fin (Paths The proof of this theorem can be found in the appendix. It remains to define the substitution of subsets of states by their abstractions. Intuitively, a subset of states is replaced by the abstraction as in Definition 9, while the incoming transitions of the initial states of the abstraction as well as the outgoing transitions of the absorbing states of the abstraction are not changed. Definition 10. Assume a PDTMC M = (S, V , I, P ), a non-absorbing set of states K ⊆ S, the induced PDTMC MK = (S K , V K , I K , P K ) and the abstraction K K K K K MK by its abstraction MK abs = (Sabs , Vabs , Iabs , Pabs ). The substitution of M abs in M is given by MK7→abs = (SK7→abs , VK7→abs , IK7→abs , PK7→abs ) with SK7→abs = K (S \ K) ∪ Sabs , VK7→abs = V and for all s, s′ ∈ SK7→abs , IK7→abs (s) = I(s) and  ′  if s ∈ /K P (s, s ) ′ K ′ PK7→abs (s, s ) = Pabs (s, s ) if s ∈ K ∧ s′ ∈ Out(K)   0 otherwise . Due to Theorem 1, it directly follows that this substitution does not change reachability properties from input states to the absorbing states of a PDTMC. Corollary 1. Given a PDTMC M and a non-absorbing subset K ⊆ S of states, it holds for all initial states sI ∈ SI and absorbing states t ∈ T that MK7→abs M Pr M (sI , t)) = Pr fin (Paths MK7→abs (sI , t)). fin (Paths 3.2 Model Checking Parametric Markov Chains In the previous section we gave the theoretical background for our model checking algorithm. Now we will describe how to compute the abstractions efficiently. As a heuristic for forming the sets of states to be abstracted, we choose an SCC-based decomposition of the graph: In Figure 1, the dashed rectangles indicate the decomposition: SCC S1 = {1, 2, 3, 4, 6, 7, 8} and the SCSs S1.1 = {2, 3, 4}, S1.2 = {6, 7, 8}, and S1.2.1 = {7, 8}. Algorithmically, Tarjan’s algorithm [12] is used to determine the SCC structure of the graph. Afterwards, for each SCC Algorithm 1 Model Checking PDTMCs abstract(PDTMC M) begin for all non-bottom SCCs K in MS\Inp(M) do K MK abs := abstract(M ) M := MK7→abs end for K := {non-absorbing states in M} M := MK7→abs return Mabs end model check(PDTMC M = (S, V , I, P ), T ⊆ S, λ ∈ Q) begin Mabs = (Sabs , Vabs , I abs , Pabs ) := abstract(M)  P P Pabs (sI , t) ≤ λ I(sI ) · return sI ∈SI (1) (2) (3) (4) (5) (6) (7) (8) (9) t∈T end K the input states Inp(K) are removed. On the resulting decomposed graph, a new search is performed, which yields a new set of SCCs which are SCSs in the original graph. This is iterated until no SCCs remain. The subset relation forms a partial order on these sets: S1.1 ⊂ S1 and S1.2.1 ⊂ S1.2 ⊂ S1 . The smallest sets according to this partial order, S1.1 and S1.2.1 , can only loop via their input state, otherwise there would be other included SCSs. Note that the deletion of the input states is only one possible heuristic for a decomposition of the graph. The general model checking algorithm is depicted in Algorithm 1. The recursive method abstract(PDTMC M) computes the abstraction Mabs by iterating over all SCCs of the graph induced by removing the input states of M (line 1). For K each SCC K, the abstraction MK is computed abs of the induced PDTMC M by a recursive call of the method (line 2, Definitions 7,9). Afterwards, MK is substituted by its abstraction inside M (line 3, Definition 10). Finally, the abstraction Mabs is computed and returned (line 7, Definition 9). This method is called by the model checking method (line 8) which yields the abstract system Mabs , in which transitions lead only from the initial states to the absorbing states. All transitions are labeled with a rational function for the reachability probability, as in Definition 6. Then the whole reachability probability is computed by building the sum of these transitions (line 9). This is compared to the given upper probability bound λ ∈ Q returning a truth-value. Note that this can be adapted for lower or strict probability bounds as well. What remains to be explained is the computation of the abstract probabilities pM abs . We distinguish the cases where the set K has one or multiple input states. One input state. We define the set of paths Rloop going from sI to sI and the set of paths Rout going from sI to some t ∈ T without revisiting sI : Rloop = {sI s1 . . . sn sI ∈ Paths M | ∀1 ≤ i ≤ n. si ∈ / {sI } ∪ T }, Rout = {sI s1 . . . sn t ∈ Paths M (1) | t ∈ T ∧ ∀1 ≤ i ≤ n. si ∈ / {sI } ∪ T }. (2) Consider a PDTMC MK induced by K with one initial state sI and the set of K i absorbing states T = {t1 , . . . , tn }. We determine the probabilities pM abs (sI , t ) for all 1 ≤ i ≤ n. As K \ Inp(K) has no non-trivial SCSs, the set Rout of outgoing paths consists of finitely many loop-free paths. The probability is computed by the following equations for all s ∈ S K :  1, if s = ti , K i P K pM (s, t ) = K ′ M ′ i abs P (s, s ) · pabs (s , t ), otherwise. (3)  s′ ∈(succ(s)∩K)\Inp(K) These probabilities can be computed by direct or indirect methods for solving linear equation systems, see, e. g., [13, Chapters 3,4]. Note that also state elimination as in [3] can be applied here. The probabilities of the abstract PDTMC MK abs = (Sabs , Iabs , Vabs , Pabs ) as in Definition 9 can now directly be computed, while an additional constraint is added in order to avoid divisions by zero:  K P  P pM MK ′ abs (sI ,t) , if MK MK (s ,t′ ) t′ ∈T pabs (sI , t ) 6= 0 p ′ I (4) Pabs (sI , t) = t ∈T abs 0, otherwise. K 1 In case there is only one absorbing state, i. e., n = 1, we have pM abs (sI , t ) = 1. This is directly exploited without further computations. Multiple input states. Given a PDTMC MK with the set of initial states SI = K i {s1I , . . . , sm I } with I (sI ) > 0 for all 1 ≤ i ≤ m and a set of absorbing states 1 n T = {t , . . . , t }. The intuitive idea would be to maintain a copy of MK for each initial state and handle the other initial states as inner states in this copy. Then, the method as described in the previous paragraph can be used. However, this is both very time and memory consuming. Therefore, we first formulate the linear K ′ ′ equation system as in Equation (3). All variables pM abs (s, s ) with s ∈ K \ Inp(K) K ′ are eliminated from the equation system. For each of the variables pM abs (sI , s ), the equation system is then solved separately by eliminating all other variables. K M Algorithm 1 returns the rational functions Pabs (sI , t) for all t ∈ T as in Equation (4). To allow only graph-preserving evaluations of the parameters, we perform preprocessing where conditions are collected according to Definition 5 as well as the ones from Equation (4). These constraints can be evaluated by a SAT modulo theories solver which can handle non-linear arithmetic over the reals [14]. In case the solver returns an evaluation which satisfies the resulting constraint set, the reachability property is satisfied. Otherwise, the property is violated. 4 Factorization of Polynomials The procedure introduced in the previous section constructs rational functions representing reachability probabilities. We now present an optimization of the frequently used arithmetic operations of addition, multiplication and division of rational functions. During the algorithm presented in Section 3 as well as the mere state-elimination [3], the rational functions that occur rapidly grow even when canceling these functions in every step. Although this exponential blow-up cannot be prevented in general, our experiments show that optimizing the arithmetic operations leads to remarkable speed ups. The key of the optimization for the operations on rational functions is to maintain a factorization for each polynomial which occurs as numerator or e e e e denominator. A polynomial g = a1 · x11,1 · . . . · xn1,n + · · · + am · x1m,1 · . . . · xnm,n is normalized if (ej,1 , . . . , ej,n ) 6= (ek,1 , . . . , ek,n ) for all j, k ∈ {1, . . . , m} with j= 6 k and the monomials are ordered, e. g., according to the reverse lexicographical ordering. A factorization Fg = {g1e1 , . . . , gnen } of a polynomial g is a set4 of factors giei , where the bases gi are normalized and pairwise different polynomials, the Qn exponents are ei ∈ N, n = 0 if g = 0, and g = i=1 giei otherwise. For polynomials g, h and a factorization Fg = {g1e1 , . . . , gnen } of g let bases(Fg ) = {g1 , . . . , gn } and exp(h, Fg ) be ei if gi = h and 0 if h ∈ / bases(Fg ). As the bases are not required to be irreducible, factorizations are not unique. We maintain that bases / Fg for g = and exponents are non-zero, F0 = ∅, F1 = {11 }, and 1k ∈ 6 1. For Fg = {g1e1 , . . . , gnen }, this is expressed by the reduction Fgred = {11 } if n > 0 and gi = 1 or ei = 0 for all 1 ≤ i ≤ n, and Fgred = Fg \ {giei | gi = 1 ∨ ei = 0} otherwise. Instead of applying arithmetic operations on two polynomials g1 and g2 directly, we operate on their factorizations Fg1 and Fg2 . We use the following operations on factorizations: Fg1 ∪F Fg2 factorizes a (not necessarily least) common multiple of g1 and g2 , Fg1 ∩F Fg2 a (not necessarily greatest) common divisor, whereas the binary operations ·F , /F , +F correspond to multiplication, division5 and addition, respectively. Due to space limitations, we omit in the remaining of this paper the trivial cases involving F0 . Therefore we define Fg1 Fg1 Fg1 Fg1 Fg1 ∪F Fg2 ∩F Fg2 ·F Fg2 /F Fg2 +F Fg2 = {hmax(exp(h,Fg1 ),exp(h,Fg2 )) | h ∈ bases(Fg1 ) ∪ bases(Fg2 )}red = {hmin(exp(h,Fg1 ),exp(h,Fg2 )) | h=1 ∨ h∈bases(Fg1 )∩bases(Fg2 )}red = {hexp(h,Fg1 )+exp(h,Fg2 ) | h ∈ bases(Fg1 ) ∪ bases(Fg2 )}red = {hmax(0,e−exp(h,Fg2 )) | he ∈ Fg1 }red = D ·F {(Πg1′ ∈Fg1 /F D g1′ ) + (Πg2′ ∈Fg2 /F D g2′ )}red where D = Fg1 ∩F Fg2 and max(a, b) (min(a, b)) equals a if a ≥ b (a ≤ b) and b otherwise. Example 3 illustrates the application of the above operations. 4 5 We represent a factorization of a polynomial as a set; however, in the implementation we use a more efficient binary search tree instead. Fg1 /F Fg2 is a factorization of g1 /g2 only if Fg1 and Fg2 are sufficiently refined and g2 divides g1 . Algorithm 2 gcd computation with factorization refinement gcd(factorization Fg1 , factorization Fg2 ) begin G := (Fg1 ∩F Fg2 ) Fi := Fgi /F G and Fi′ := {11 } for i = 1, 2 while exists r1e1 ∈ F1 with r1 6= 1 do F1 := F1 \ {r1e1 } while r1 6= 1 and exists r2e2 ∈ F2 with r2 6= 1 do F2 := F2 \ {r2e2 } if ¬irreducible(r1 ) ∨ ¬irreducible(r2 ) then g := common gcd(r1 , r2 ) else g := 1 if g = 1 then F2′ := F2′ ·F {r2e2 } else r1 := rg1 Fi := Fi ·F {g ei −min(e1 ,e2 ) } for i = 1, 2 F2′ := F2′ ·F {( rg2 )e2 } G := G ·F {g min(e1 ,e2 ) } end if end while F1′ := F1′ ·F {r1e1 } F2 := Fg2′ F2′ := {11 } end while return (F1′ , F2 , G) end (1) (2) (3) (4) (5) (6) (7) (8) (9) (10) (11) (12) (13) (14) (15) (16) (17) (18) (19) (20) (21) (22) For rational functions hg we maintain separate factorizations Fg and Fh for the numerator g and the denominator h, respectively. For multiplication hg = hg11 · hg22 , we compute Fg = Fg1 ·F Fg2 and Fh = Fh1 ·F Fh2 . For division we use the multiplication due to hg11 : hg22 = hg11 · hg22 . For the addition hg = hg11 + hg22 , we compute h as a common multiple of h1 and h2 factorized by Fh = Fh1 ∪F Fh2 , such that h = hi · h′i with Fh′i = Fh /F Fhi for i = 1, 2. For the numerator g we first determine a common divisor d of g1 and g2 by Fd = Fg1 ∩F Fg2 , such that gi = d · gi′ with Fgi′ = Fgi /F Fd for i = 1, 2. The numerator g is d·(g1′ ·h′1 +g2′ ·h′2 ) with factorization Fd ·F (Fg1′ ·F Fh′1 +F Fg2′ ·F Fh′2 ). The rational function hg resulting from the addition is further simplified by cancelation, i. e., dividing g and h by their greatest common divisor (gcd) g ′ . Given the factorizations Fg and Fh , Algorithm 2 calculates the factorizations Fg′ , F g1′ , and F g2′ by reusing Fg and Fh as much as possible. Initially, a factorization g g G of a common divisor of g1 and g2 is set to Fg1 ∩F Fg2 . The factors of g1 and g2 , which are not part of the factorization, are stored in F1 resp. F2 . Moreover, F1′ and F2′ contain those factors of g1 resp. g2 , such that for all (f1′ , f2′ ) ∈ F1′ × F2′ f1′ and f2′ have no common divisors. The algorithm now iteratively adds further common divisors of g1 and g2 to G until it is a factorization of their gcd. For this purpose, we consider for each factor in F1 all factors in F2 and calculate the gcd of their bases which we compute by a common approach for calculating gcds. Note that the main concern of Algorithm 2 is to avoid this rather expensive operation or else call it on preferably simple polynomials. Where the latter is entailed by the idea of using factorizations, the former can be achieved by excluding pairs of factors for which we can cheaply decide that both are irreducible, i. e., they have no non-trivial divisors. If factors r1e1 ∈ F1 and r2e2 ∈ F2 with g := common gcd(r1 , r2 ) = 1 are found, we just shift r2e2 from F2 to F2′ . Otherwise, we can add g min(e1 ,e2 ) , which is the gcd of r1e1 and r2e2 , to G and extend the factors F1 resp. F2 , which could still contain common divisors, by g e1 −min(e1 ,e2 ) resp. g e2 −min(e1 ,e2 ) . Furthermore, F2′ obtains the new factor ( rg2 )e2 , which has certainly no common divisor with any factor in F1′ . Finally, we set the basis r1 to rg1 , excluding the just found common divisor. If all factors in F2 have been considered for common divisors with r1 , we can add it to F1′ and continue with the next factor in F1 , for which we must reconsider all factors in F2′ and, therefore, shift them to F2 . The algorithm terminates, if the last factor of F1 has been processed, returning the factorizations Fg′ , F g1′ and F g2′ , which we can use to refine the factorizations of g1 and g2 via g g Fg1 := F g1′ ·F G and Fg2 := F g2′ ·F G. g g Example 3. Assume we want to apply Algorithm 2 to the factorizations Fxyz = {(xyz)1 } and Fxy = {(x)1 , (y)1 }. We initialize G = F1′ = F2′ = {(1)1 }, F1 = Fxyz and F2 = Fxy . First, we choose the factors (r1 )e1 = (xyz)1 and (x)1 and remove them from F1 resp. F2 . The gcd of their bases is x, hence we only update r1 to (yz)1 and G to {(x)1 }. Then we remove the next and last element (y)1 from F2 . Its basis and r1 have the gcd y and we therefore update r1 to (z)1 and G to {(x)1 , (y)1 }. Finally, we add (z)1 to F1′ and return the expected result ({(z)1 }, {(1)1 }, {(x)1 , (y)1 }). Furthermore, we can update Fxyz = F1′ ·F G = {(x)1 , (y)1 , (z)1 } afterwards. Theorem 2. Let p1 and p2 be polynomials with factorizations Fp1 resp. Fp2 as before. Applying Algorithm 2 to these factorizations results in gcd(Fp1 , Fp2 ) = (Fr1 , Fr2 , G) with G being a factorization of the greatest common divisor g of p1 and p2 and Fr1 and Fr2 being factorizations of pg1 resp. pg2 . The proof of this theorem can be found in the appendix. 5 Experiments We developed a C++ prototype implementation of our approach using the arithmetic library GiNaC [15]. The prototype is available for testing on the project homepage6 . Moreover, we implemented the state-elimination approach that is used by PARAM [6] using our optimized factorization approach to provide a more distinct comparison. All experiments were run on an Intel Core 2 Quad CPU 2.66 GHz with 4 GB of memory. We defined a timeout (T O) of 3600 seconds and a 6 http://goo.gl/nS378q memory bound (M O) of 4 GB. We report on three case studies; a more distinct description and the specific instances we used are available at our homepage. The bounded retransmission protocol (BRP) [16] models the sending of files via an unreliable network, manifested in two lossy channels for sending and acknowledging the reception. This model is parametrized in the probability of reliability of those channels. The crowds protocol (CROWDS) [17] is designed for anonymous network communication using random routing, parametrized in how many members are “good” or “bad” and if a good member delivers a message or randomly routes it to another member. NAND multiplexing (NAND) [18] models how reliable computations are obtained using unreliable hardware by having a certain number of copies of a NAND unit all doing the same job. Parameters are the probabilities of faultiness of the units and of erroneous inputs. The experimental setting includes our SCC-based approach as described in Section 3 using the optimized factorization of polynomials as in Section 4 (SCC MC), the state elimination as in PARAM but also using the approach of Section 4 (STATE ELIM) and PARAM itself.7 For all instances we list the number of states and transitions; for each tool we give the running time in seconds and the memory consumption in MB; the best time is boldfaced. Moreover, for our approaches we mention the number of polynomials which are intermediately stored. Model — Graph States Trans. SCC MC Time Poly STATE ELIM Mem Time Poly PARAM Mem Time Mem BRP 2695 3459 1.86 2319 16.04 1.97 6647 40.04 36.49 17.38 BRP 5383 6915 14.81 4623 47.00 12.19 13367 161.80 356.04 63.21 BRP 10378 13827 TO — — 63.47 21805 328.94 3203.20 431.40 BRP 10759 13827 147.31 9231 176.89 85.54 26807 682.24 3511.96 304.07 CROWDS 8655 14953 4.15 8747 13.21 3.24 2943 11.96 139.70 10.44 CROWDS 37293 65011 16.69 33549 40.23 21.72 8148 30.61 1977.95 35.39 CROWDS 198201 348349 80.05 155000 175.40 271.04 27344 133.82 TO — CROWDS 726379 1283297 262.88 500048 668.94 1742.42 73702 477.26 TO — CROWDS 2888763 5127151 1030.57 1707776 2705.35 TO — — TO — NAND 14323 21567 39.71 25504 366.79 59.60 405069 926.33 15.26 16.89 NAND 28183 42287 208.41 44799 1405.16 218.85 925324 3708.27 50.45 30.47 NAND 35113 52647 352.09 54445 2047.66 364.09 1184848 3696.39 78.19 40.51 NAND 78334 121512 639.29 184799 3785.11 TO — — 1127.83 113.32 For BRP, STATE ELIM outperforms PARAM by up to two orders of magnitude while SCC MC runs into a timeout for one particular instance. This is due to the graph structure. In contrast, the crowds protocol induces a nested SCC structure, which is very hard for PARAM since many divisions of polynomials have to be carried out. On larger benchmarks, it is therefore outperformed by more than three orders of magnitude while SCC MC performs best. The NAND graphs consist of single paths, which induces a high number of polynomials we are keeping for the factorization. Our implementation offers the possibility to 7 Note that no bisimulation reduction was applied to any of the input models, which would improve the feasibility of all approaches likewise. bound this pool of polynomials which highly decreases the memory consumption for the sake of loosing information about the factorizations. 6 Conclusion and Future Work We presented a new approach to verify parametric Markov chains together with an improved factorization of polynomials. We were able to highly improve the scalability in comparison to existing approaches. Future work will be dedicated to the actual parameter synthesis. First, we want to incorporate interval constraint propagation [19] in order to provide reasonable intervals for the parameters where properties are satisfied or violated. Moreover, we are going to investigate the possibility of extending our approaches to models with costs. References 1. Su, G., Rosenblum, D.S.: Asymptotic bounds for quantitative verification of perturbed probabilistic systems. In: Proc. of ICFEM. Volume 8144 of LNCS, Springer (2013) 297–312 2. Daws, C.: Symbolic and parametric model checking of discrete-time Markov chains. In: Proc. of ICTAC. Volume 3407 of LNCS, Springer (2004) 280–294 3. Hahn, E.M., Hermanns, H., Zhang, L.: Probabilistic reachability for parametric Markov models. Software Tools for Technology Transfer 13(1) (2010) 3–19 4. Gruber, H., Johannsen, J.: Optimal lower bounds on regular expression size using communication complexity. In: Proc. of FOSSACS. Volume 4962 of LNCS, Springer (2008) 273–286 5. Hahn, E.M., Han, T., Zhang, L.: Synthesis for PCTL in parametric Markov decision processes. In: Proc. of NFM. Volume 6617 of LNCS, Springer (2011) 146–161 6. Hahn, E.M., Hermanns, H., Wachter, B., Zhang, L.: PARAM: A model checker for parametric markov models. In: Proc. of CAV. Volume 6174 of LNCS, Springer (2010) 660–664 7. Han, T., Katoen, J.P., Mereacre, A.: Approximate parameter synthesis for probabilistic time-bounded reachability. In: Proc. of RTSS, IEEE CS (2008) 173–182 8. Puggelli, A., Li, W., Sangiovanni-Vincentelli, A.L., Seshia, S.A.: Polynomial-time verification of PCTL properties of MDPs with convex uncertainties. In: Proc. of CAV. Volume 8044 of LNCS, Springer (2013) 527–542 9. Benedikt, M., Lenhardt, R., Worrell, J.: LTL model checking of interval Markov chains. In: Proc. of TACAS. Volume 7795 of LNCS, Springer (2013) 32–46 10. Baier, C., Katoen, J.P.: Principles of Model Checking. The MIT Press (2008) 11. Ábrahám, E., Jansen, N., Wimmer, R., Katoen, J.P., Becker, B.: DTMC model checking by SCC reduction. In: Proc. of QEST, IEEE CS (2010) 37–46 12. Tarjan, R.E.: Depth-first search and linear graph algorithms. SIAM Journal on Computing 1(2) (1972) 146–160 13. Quarteroni, A., Sacco, R., Saleri, F.: Numerical Mathematics. Springer (2000) 14. Jovanovic, D., de Moura, L.M.: Solving non-linear arithmetic. In: Proc. of IJCAR. Volume 7364 of LNCS, Springer (2012) 339–354 15. Bauer, C., Frink, A., Kreckel, R.: Introduction to the GiNaC framework for symbolic computation within the C++ programming language. J. Symb. Comput. 33(1) (2002) 1–12 16. Helmink, L., Sellink, M., Vaandrager, F.: Proof-checking a data link protocol. In: Proc. of TYPES. Volume 806 of LNCS, Springer (1994) 127–165 17. Reiter, M.K., Rubin, A.D.: Crowds: Anonymity for web transactions. ACM Trans. on Information and System Security 1(1) (1998) 66–92 18. Han, J., Jonker, P.: A system architecture solution for unreliable nanoelectronic devices. IEEE Transactions on Nanotechnology 1 (2002) 201–208 19. Fränzle, M., Herde, C., Teige, T., Ratschan, S., Schubert, T.: Efficient solving of large non-linear arithmetic constraint systems with complex boolean structure. Journal on Satisfiability, Boolean Modeling, and Computation 1(3-4) (2007) 209–236 Appendix Here, we show the correctness of the PDTMC abstraction as in Definition 9. First, let us recall Theorem 1 Theorem 3. For a PDTMC M = (S, V , I, P ) and its abstraction Mabs = (Sabs , Iabs , Vabs , Pabs ) according to Definition 9 it holds for all initial states sI ∈ SI and all absorbing states t ∈ T that M Mabs abs Pr M (sI , t)) = Pr M (sI , t)). fin (Paths fin (Paths Proof. As the bottom SCCs are exactly the absorbing target states in T , the probability of reaching a state of T is 1. The probability pM abs (sI , sI ) can therefore be expressed w. r. t. the probabilities of reaching an absorbing state without revisiting sI : X pM (5) pM abs (sI , t). abs (sI , sI ) = 1 − t∈T To reduce notation, we define the set of paths Rloop looping on sI and the set of paths Rout going to some t ∈ T without revisiting sI . Rloop ={sI s1 . . . sn sI ∈ Paths M | si ∈ / {sI } ∪ T, 1 ≤ i ≤ n} (6) M (7) Rout ={sI s1 . . . sn t ∈ Paths | si ∈ / {sI } ∪ T, 1 ≤ i ≤ n, t ∈ T } As the self-loop in sI represents the paths of Rloop , it holds that pM abs (sI , sI ) = Pr fin (Rloop ). We now have: = M Pr M (sI , t)) fin (Paths ∞ [ {π1 · · · · · πi · πout | πj ∈ Rloop , 1 ≤ j ≤ i; πout ∈ Rout }) Pr M fin ( i=0 = ∞ X Pr M fin ({π1 · · · · · πi · πout | πj ∈ Rloop , 1 ≤ j ≤ i; πout ∈ Rout }) = ∞ X M i (Pr M fin (Rloop )) · Pr fin (Rout ) = ∞ X M i (pM abs (sI , sI )) · Pr fin (Rout ) (Equation (8)) i=0 i=0 i=0 = = 1 1− pM abs (sI , sI ) · Pr M fin (Rout ) (Geometric Series) 1 P sout ∈T pM abs (sI , sout ) · Pr M fin (Rout ) (Equation (5)) (8) = 1 P sout ∈T pM abs (sI , sout ) · pM abs (sI , t) (Definition 8) = Pabs (sI , t) (Definition 9) = Mabs abs Pr M (sI , t)) fin (Paths As the probabilities of reaching the absorbing states from initial states coincide in M and Mabs , our abstraction is valid. We show the correctness and completeness of Algorithm 2 by proving the following theorem. Theorem 4. Let p1 and p2 be polynomials with factorizations Fp1 resp. Fp2 as before. Applying Algorithm 2 to these factorizations results in gcd(Fp1 , Fp2 ) = (Fr1 , Fr2 , G) with G being a factorization of the greatest common divisor g of p1 and p2 and Fr1 and Fr2 being factorizations of pg1 resp. pg2 . Proof. We define the product of a factorization Fp by P(Fp ) = the common gcd by gcdc . We first prove that F1 ·F F1′ ·F {(r1 )e1 } ·F G = Fp1 , F2 ·F F2′ ·F G = Fp2 , gcdc (P(F1′ ·F {(r1 )e1 }), P(F2′ )) = 1 Q q e ∈Fp q e and (9) (10) (11) hold after the ith pass through the inner while-loop (line 16) of Algorithm 2, for all i ∈ N ∪ {0}, where we can assume without loss of generality that (r1 )e1 is initialized by (1)1 in the beginning of Algorithm 2. In the following we denote by (i) ′(i) (i) ′(i) (i) F1 , F1 , F2 , F2 , r1 , G(i) the according datastructures in the ith pass through the inner while-loop. Basis (i = 0): Before entering the outer while-loop and, hence, before entering the inner while-loop (line 2) Equation 9 and Equation 10 hold as a consequence ′(0) ′(0) of the definiton of ·F . Equation 11 is trivially implied from F1 = F2 = 1 {(1) }. Inductive step (i > 0): If g = 1 (line 9) we only shift r2e2 from F2 to F2′ , thus the left-hand sides of Equation 9 and Equation 10 have the same result as in the previous pass through the inner while-loop (i−1) and both equations follow from the inductive hypothesis. Concerning Equation 11 either F1′ = {(1)1 } or all elements it contains have been added directly after exiting the inner while-loop and therefore the inductive hypothesis still holds. After resetting F2 in line 19 it holds that gcdc (P(F1′ ·F {(r1 )e1 }), P(F2 )) = gcdc (P(F1′ ·F {(r1 )e1 }), P(F2′ )) IH = 1 As we do not add elements to F2 nor F1′ anywhere else, P(F2 ) and P(F1′ ) are always coprime (*), which means that they have no common divisors. It follows that ′(i) gcdc (P(F1 ′(i) ·F {(r1,i )e1 }), P(F2 )) ′(i−1) ′(i−1) ·F {(r1,i )e1 }), P(F2 = gcdc (P(F1 ′(i−1) ) · (r1,i )e1 , P(F2 ′(i−1) ) · (r1,i )e1 , P(F2 ′(i−1) ) · (r1,i )e1 , r2e2 ) ′(i−1) ), r2e2 ) · gcdc ((r1,i )e1 , r2e2 ) = gcdc (P(F1 ≤ gcdc (P(F1 IH ·F {(r2 )e2 })) = gcdc (P(F1 ≤ gcdc (P(F1 ′(i−1) ) · r2e2 ) ′(i−1) )) · gcdc (P(F1 ′(i−1) ) · (r1,i )e1 , r2e2 ) (∗) = gcdc ((r1,i )e1 , r2e2 ) g=1 = 1 and therefore Equation 11 holds. If g 6= 1, then Equation 9 holds as (i) ′(i) F1 · F F1 (i−1) = F1 ·F {(r1,i )e1 } ·F G(i) ′(i−1) ·F {g e1 −min(e1 ,e2 ) } ·F F1 ·F {( r1,i−1 e1 ) } g ·F G(i−1) ·F {g min(e1 ,e2 ) } (i−1) = F1 ′(i−1) · F F1 ·F {(r1,i−1 )e1 } ·F G(i−1) IH = Fp1 and Equation 10 holds because of (i) ′(i) F2 · F F2 (i−1) = (F2 ·F G(i) ′(i−1) /F {r2e2 }) ·F {g e2 −min(e1 ,e2 ) } ·F F2 ·F G(i−1) ·F {g min(e1 ,e2 ) } (i−1) = F2 IH = Fp2 . ′(i−1) · F F2 ·F G(i−1) ·F {( r2 e2 ) } g Furthermore, Equation 11 holds as a consequence of ′(i) ′(i) ·F {(r1,i )e1 }), P(F2 )) r1,i−1 e1 r2 ′(i−1) ′(i−1) = gcdc ( P(F1 ·F {( ) }), P(F2 ·F {( )e2 }) ) g g gcdc (P(F1 ′(i) ′(i−1) ·F {(r1,i−1 )e1 }), P(F2 )) r r2 r 1,i−1 e1 2 ′(i−1) ) , ( ) e2 ) · gcdc (P(F1 ), ( )e2 ) · gcdc (( g g g r1,i−1 e1 r2 e2 r 2 e2 IH ′(i−1) ) , ( ) ) = gcdc (P(F1 ), ( ) ) · gcdc (( g g g r 2 e2 g gcd ′(i−1) ), ( ) ) = gcdc (P(F1 of r1 r2 g ≤ gcdc (P(F1 (∗) = 1. Now we can prove the completeness and correctness of Algorithm 2. Completeness: When passing the outer while-loop it holds that ′(i) P(F1 ) ′(i−1) P(F1 ) ′(i) ≥ 1, P(F2 ) ′(i−1) P(F2 ) ≥ 1 and P(G(i) ) ≥ 1, P(G(i−1) ) and for at least one of these inequalities the relation is even strict (>). From Equation 9 and Equation 10, the fact that Fp1 and Fp2 are constant and we only consider factorizations of polynomials ≥ 1 (in particular F1 and F2 ) follows the completeness of Algorithm 2. Correctness: When Algorithm 2 leaves the outer while-loop it holds that F1 = F2′ = {(1)1 }. This and the fact that Equation 9 and Equation 10 are valid in the end of the inner while-loop implies: F1′ ·F G = Fp1 F2 ·F G = Fp2 As also Equation 11 holds, Algorithm 2 is correct.