Hongming Liu
email: hm-liu@sjtu.edu.cn

Shanghai Jiao Tong

Hongfei Fu

Zhiyong Yu
email: yuzhiyong18@sjtu.edu.cn

Jiaxin Song

Guoqiang Li
email: li.g@sjtu.edu.cn

Scalable Linear Invariant Generation with Farkas' Lemma

Invariant generation is the classical problem to automatically generate invariants to aid the formal analysis of programs. In this work, we consider the linear-invariant-generation problem over affine programs (i.e., programs with affine guards and updates). In the literature, the only known sound and complete characterization to solve this problem is via Farkas' Lemma (FL), and has been implemented through either quantifier elimination or reasonable heuristics. Although FL-based approaches can generate highly accurate linear invariants from the completeness of FL, the main bottleneck to apply these approaches is the scalability issue caused by either non-linear constraints or combinatorial explosion. It has long been an unresolved problem to improve the scalability of FL-based approaches. In this work, we address this issue with novel algorithmic improvements. In detail, we base our approach on the only practical FL-based approach [Sankaranarayanan et al., SAS 2004] that applies FL with reasonable heuristics, and develop two novel and independent improvements to leverage the scalability. The first improvement is the novel idea to generate invariants at only one program location in a single invariant-generation process, so that the invariants for each location are generated separately rather than together in a single computation. This idea naturally leads to a parallel processing that divides the invariant-generation task for all program locations by assigning the locations separately to mulitple processors. Moreover, the idea enables us to develop detailed technical improvements to further reduce the combinatorial explosion in the original work [Sankaranarayanan et al., SAS 2004]. The second improvement is a segmented subsumption testing in the CNF-to-DNF expansion that allows to discover more local subsumptions in advance. We formally prove that our approach has the same accuracy as the original work, thus does not incur accuracy loss on the generated invariants. Moreover, experimental results on representative benchmarks involving non-trivial linear invariants demonstrate that our approach improves the runtime of the original work by several orders of magnitude, even if in the non-parallel scenario that sums up the execution time for all program locations. Thus, our approach constitutes the first significant improvement in FL-based approaches for linear invariant generation after almost two decades.

INTRODUCTION

Invariants. An assertion at a program location is called an invariant if it is always satisfied by the values taken by the program variables whenever the location is reached in the execution of the program. Invariants play a fundamental role in program analysis and verification as they provide over-approximation for reachable program states. Therefore, they are widely used in proving basic properties such as safety [START_REF] Albarghouthi | Ufo: A Framework for Abstraction-and Interpolation-Based Software Verification[END_REF][START_REF] Manna | Temporal verification of reactive systems -safety[END_REF][START_REF] Padon | Ivy: safety verification by interactive generalization[END_REF], reachability [START_REF] Alias | Multi-dimensional Rankings, Program Termination, and Complexity Bounds of Flowchart Programs[END_REF][START_REF] Asadi | Polynomial reachability witnesses via Stellensätze[END_REF][START_REF] Bradley | Linear Ranking with Reachability[END_REF][START_REF] Chen | Discovering Non-linear Ranking Functions by Solving Semi-algebraic Systems[END_REF][START_REF] Colón | Synthesis of Linear Ranking Functions[END_REF][START_REF] David | Danger Invariants[END_REF][START_REF] Podelski | A Complete Method for the Synthesis of Linear Ranking Functions[END_REF] and time-complexity analysis [START_REF] Chatterjee | Non-polynomial Worst-Case Analysis of Recursive Programs[END_REF], etc. The quality of the generated invariants is measured by their accuracy, i.e., the amount of over-approximation against the actual set of reachable program states. The accuracy of the invariants is an important factor as inaccurate invariants can lead to loose results or even failure to get meaningful results for program analysis and verification. Invariant Generation. Invariant generation is the classical problem that asks to automatically generate invariants for an input program, and has been studied for decades. Various approaches have been proposed to solve the problem, such as abstract interpretation [START_REF] Cousot | Abstract Interpretation: A Unified Lattice Model for Static Analysis of Programs by Construction or Approximation of Fixpoints[END_REF][START_REF] Cousot | Automatic Discovery of Linear Restraints Among Variables of a Program[END_REF], constraint solving [START_REF] Chatterjee | Polynomial invariant generation for non-deterministic recursive programs[END_REF][START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Kapur | Automatically Generating Loop Invariants Using Quantifier Elimination[END_REF], recurrence analysis [START_REF] Farzan | Compositional Recurrence Analysis[END_REF][START_REF] Humenberger | Algebra-Based Synthesis of Loops and Their Invariants[END_REF][START_REF] Kincaid | Compositional recurrence analysis revisited[END_REF][START_REF] Kincaid | Non-linear reasoning for invariant synthesis[END_REF], logical inference [START_REF] Dillig | Inductive invariant generation via abductive inference[END_REF][START_REF] Gan | Nonlinear Craig Interpolant Generation[END_REF][START_REF] Garg | ICE: A Robust Framework for Learning Invariants[END_REF][START_REF] Gulwani | Constraint-Based Invariant Inference over Predicate Abstraction[END_REF][START_REF] Kenneth | Quantified Invariant Generation Using an Interpolating Saturation Prover[END_REF][START_REF] Sharma | From invariant checking to invariant inference using randomized search[END_REF], machine learning [START_REF] Garg | Learning invariants using decision trees and implication counterexamples[END_REF][START_REF] He | Learning fast and precise numerical analysis[END_REF][START_REF] Yao | Learning nonlinear loop invariants with gated continuous logic networks[END_REF], dynamic analysis [START_REF] Csallner | DySy: dynamic symbolic execution for invariant inference[END_REF][START_REF] Nguyen | Using dynamic analysis to discover polynomial and array invariants[END_REF][START_REF] Sharma | A Data Driven Approach for Algebraic Loop Invariants[END_REF], etc.

To infer invariants at program locations directly is often infeasible, and most existing approaches generate invariants by considering a strengthened notion called inductive invariants. An inductive invariant at a program location is an assertion that holds for the first visit to the location and is preserved under every cyclic execution path to and from the location. Inductive invariants are guaranteed to be invariants, and the well-established method to prove that an assertion is an invariant is to find an inductive invariant that strengthens it [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Manna | Temporal verification of reactive systems -safety[END_REF].

Numerical Inductive Invariants. An important category of inductive invariants is that of numerical inductive invariants that captures the relationship between the numerical values taken by the program variables. Numerical values are a basic aspect of programs, and many common failures of programs (such as array out-of-bound, division by zero, etc.) are closely related to numerical values. Thus, numerical inductive invariants are essential in proving numeric-critical properties. In this work, we consider automated generation of numerical inductive invariants. To be more precise, we consider algorithmic approaches for generating linear inductive invariants below.

A notable subclass of numerical inductive invariants is the class of linear inductive invariants. Informally, a numerical inductive invariant is linear if the invariant takes the form of a system of linear inequalities over the program variables. Linear inductive invariants are the most basic form of numerical invariant invariants, hence is important both for the academic and practical purpose. To resolve the automated generation of linear inductive invariants, we consider the method of constraint solving, as follows.

The Method of Constraint Solving. To solve the invariant-generation problem, constraint-solving based approaches usually consider the following paradigm: first establish a template with unknown parameters for the target invariant, then collect constraints from the inductive condition for invariants, and finally solve the unknown parameters in the template to get the desired invariants.

Constraint-solving based approaches for numerical invariant generation can roughly be classified by linear and polynomial invariant generation. For linear invariant generation, Farkas' Lemma provides a complete characterization of the inductive condition and has been studied in [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], which was further solved by quantifier elimination [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF] and several heuristics [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. The STING invariant generator [START_REF] Sting | StInG: Stanford Invariant Generator[END_REF] implements the approach in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], and the INVGEN invariant generator [START_REF] Gupta | InvGen: An Efficient Invariant Generator[END_REF] integrates abstraction interpretation and the approach in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. Besides, an approach based on eigenvectors for a restricted class of invariants is proposed in [START_REF] Steven De Oliveira | Synthesizing Invariants by Solving Solvable Loops[END_REF]. Recently, probabilistic linear invariants have also been considered in probabilistic programs through Farkas' Lemma and Motzkin's Transposition theorem [START_REF] Chatterjee | Stochastic invariants for probabilistic termination[END_REF][START_REF] Katoen | Linear-Invariant Generation for Probabilistic Programs: -Automated Support for Proof-Based Methods[END_REF].

For polynomial invariant generation, a variety of approaches were proposed in the literature. Complete approaches (that typically have high runtime complexity) were proposed through quantifier elimination [START_REF] Kapur | Automatically Generating Loop Invariants Using Quantifier Elimination[END_REF] and other computer-algebra based techniques [START_REF] Yang | Recent advances in program verification through computer algebra[END_REF]. Semi-complete approaches (that have lower runtime complexity but retain completeness in restricted situations as compared with complete approaches) through Positivstellensätze have been proposed in [START_REF] Chatterjee | Polynomial invariant generation for non-deterministic recursive programs[END_REF]. The special case of polynomial-equality invariants was solved completely by Zariski-closure [START_REF] Hrushovski | Polynomial Invariants for Affine Programs[END_REF] and Gröbner basis [START_REF] Rodríguez | Automatic Generation of Polynomial Loop Invariants: Algebraic Foundations[END_REF]. Heuristics for polynomial invariant generation have also been extensively studied, such as semidefinite programming with relaxation [START_REF] Adjé | Property-based Polynomial Invariant Generation Using Sums-of-Squares Optimization[END_REF][START_REF] Cousot | Proving Program Invariance and Termination by Parametric Abstraction, Lagrangian Relaxation and Semidefinite Programming[END_REF][START_REF] Lin | Proving total correctness and generating preconditions for loop programs via symbolic-numeric computation methods[END_REF], Lagrange interpolation [START_REF] Chen | Counterexample-Guided Polynomial Loop Invariant Generation by Lagrange Interpolation[END_REF], reduction to linear algebra [START_REF] Steven De Oliveira | Polynomial Invariants by Linear Algebra[END_REF], Hypergeometric sequences [START_REF] Humenberger | Automated Generation of Non-Linear Loop Invariants Utilizing Hypergeometric Sequences[END_REF] and Gröbner basis [START_REF] Sriram Sankaranarayanan | Non-linear loop invariant generation using Gröbner bases[END_REF]. Recently, an approach based on Stengle's Postivstellensatz for generating probabilistic polynomial invariants in probabilistic programs is proposed in [START_REF] Feng | Finding Polynomial Loop Invariants for Probabilistic Programs[END_REF].

Compared with other methods (such as abstract interpretation, machine learning, etc.), constraint solving has the advantage of a theoretical guarantee on the accuracy of the generated invariants based on the considered numerical form of the invariants and the accuracy loss caused by the heuristics (if any), but typically require higher runtime complexity.

Our Contribution. We follow the constraint-solving method and consider automated generation of linear inductive invariants over affine programs. An affine program is an imperative program where every assignment is an affine expression over program variables and every guard condition (in e.g. conditional branches, while loops, etc.) is a propositional combination of comparison between affine expressions over program variables. Our work is based on the previous work [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that implements FL with reasonable heuristics and is the only FL-based approach that achieves practical performance over realistic programs.

We propose two novel improvements to leverage the scalability of applying FL to linear invariant generation, leading to a novel approach that significantly improves the scalability of [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. In [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], the scalability issue arises from two sources of combinatorial explosion: a CNF-to-DNF expansion and also a polyhedral-cone generator computation. Our first improvement is the novel idea that instead of solving the invariants at all program locations in one invariant-generation process when applying FL, our approach solves the invariants at each location separately, i.e., solving the invariants location-by-location. The idea directly leads to the potential to speed up the invariantgeneration process at the target location since one does not need to generate the invariants at other locations, and naturally enables further speed up through parallel processing (that computes the invariants for every location separately over multiple processors). Moreover, the idea enables us to develop the following technical improvements to further mitigate the combinatorial explosion in the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]: (i) in the CNF-to-DNF expansion, we reorder the expansion so that the target location (over which the invariants are to be generated) comes first to detect subsumptions earlier; (ii) in the polyhedral-cone generator computation, we eliminate variables unrelated to the target location in advance to speed up the generator computation. Our second improvement is the extra segmented subsumptiong testing that serves as a pre-processing for the CNF-to-DNF expansion and could detect local subsumption in advance.

We develop our approach on the abstract model of linear transitions systems [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that capture general affine updates (with affine guards) between program locations, so that our approach is applicable to general affine programs. To complement the generality of our approach, we implement a prototype transformation from C program codes into linear transition systems to showcase the practical connection with concrete programs. We show that our approach generates exactly the same linear invariants as the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] does, thus inherits the merit of high accuracy from constraint solving. Moreover, experimental results on typical benchmarks which involves non-trivial linear invariants reveal that our approach can indeed substantially improves the scalability against the state-of-the-art approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. For example, for several large-scale benchmarks our approach is able to attain up to thousands of orders of magnitude speed-up when compared with [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. As a result, our approach constitutes to our best knowledge the first significant improvement in Farkas'-Lemma based methods after almost two decades (since [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]).

LINEAR TRANSITION SYSTEMS AND INVARIANTS

We consider linear transition systems (LinTS's) [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] as the underlying model for invariant generation. A LinTS is composed of locations and linear transitions between locations, thus is suitable for modelling the executions of an affine program, for which a location in a LinTS corresponds to a program location (a.k.a program counter) of an affine program, and the linear transitions corresponds to the affine updates (arising from assignment statements) and guards (from if-branches and while-loops) in the program.

To present the definitions for LinTS's, we first define the basic notions of linear (in)equalities and assertions. For linear inequalities, we only consider the non-strict comparison operator ≥. Note that although an equality 𝛼 = 𝛽 can be equivalently expressed by two inequalities 𝛼 ≤ 𝛽 and 𝛼 ≥ 𝛽, the equalities in a LinTS are tackled directly as various optimizations could be applied to equalities. Also note that inequalities of the form 𝛼 ≤ 𝛽 could be equivalently tranformed into -𝛼 ≥ -𝛽. Linear (In)equalities and Assertions. A linear equality over a set 𝑉 = {𝑥 1 , . . . , 𝑥 𝑛 } of real-valued variables is of the form 𝑎 1 𝑥 1 + • • • + 𝑎 𝑛 𝑥 𝑛 + 𝑏 = 0, where 𝑎 𝑖 's and 𝑏 are real coefficients. Analogously, a linear inequality over 𝑉 is of the form 𝑎 1 𝑥 1 + • • • + 𝑎 𝑛 𝑥 𝑛 + 𝑏 ≥ 0. A linear assertion over 𝑉 is a conjunction of linear equalities and inequalities over 𝑉 .

Then we present the abstract model of linear transition systems. Linear Transition Systems. A linear transition system (LinTS) is a tuple ⟨𝑋, 𝑋 ′ , 𝐿, T, ℓ * , 𝜃 ⟩ where:

• 𝑋 is a finite set of real-valued variables such that each variable 𝑥 ∈ 𝑋 represents the current value of the variable, while 𝑋 ′ = {𝑥 ′ | 𝑥 ∈ 𝑋 } is the corresponding set of primed variables such that each primed variable 𝑥 ′ ∈ 𝑋 ′ represents the value of the unprimed counterpart 𝑥 ∈ 𝑋 in the next step of the system; • 𝐿 is a finite set of locations and ℓ * ∈ 𝐿 is the initial location; • T is a finite set of transitions where each transition 𝜏 is a triple ⟨ℓ, ℓ ′ , 𝜌⟩ that specifies the jump from the current location ℓ to the next location ℓ ′ with the guard condition 𝜌 as a linear assertion over 𝑋 ∪ 𝑋 ′ ; • 𝜃 is a linear assertion over the variables 𝑋 that specifies the initial condition at the initial location ℓ * .

To describe the behaviour of a LinTS, we further define the notions of valuations, configurations and their associated satisfaction relation as follows.

Valuations and Configurations. A valuation over a variable set 𝑉 is a function 𝜎 : 𝑉 → R that assigns to each variable 𝑥 ∈ 𝑉 a real value 𝜎 (𝑥) that corresponds to the current value held by 𝑥. In this work, we mainly consider valuations over the variable set 𝑋 of a LinTS and simply abbreviate "valuation over 𝑋 " as "valuation" (i.e., omitting 𝑋 ). Given a LinTS, a configuration is a pair (ℓ, 𝜎) such that ℓ ∈ 𝐿 is a location and 𝜎 is a valuation (over 𝑋 ), with the intuition that ℓ is the current location and 𝜎 specifies the current values for the variables in the LinTS. The Satisfaction Relation. Given a linear assertion 𝜑 over 𝑋 and a valuation 𝜎, we write 𝜎 |= 𝜑 to mean that 𝜎 satisfies 𝜑, i.e., 𝜑 is true when one substitutes the corresponding values 𝜎 (𝑥) in 𝜎 to all the variables 𝑥 in 𝜑. Analogously, given two valuations 𝜎, 𝜎 ′ (over 𝑋 ) and a linear assertion 𝜑 over 𝑋 ∪ 𝑋 ′ , we write 𝜎, 𝜎 ′ |= 𝜑 to mean that 𝜑 is true when one substitutes every variable 𝑥 ∈ 𝑋 by 𝜎 (𝑥) and every variable 𝑥 ′ ∈ 𝑋 ′ by 𝜎 ′ (𝑥) in 𝜑. Moreover, given two linear assertions 𝜑,𝜓 over 𝑋 , we write 𝜑 |= 𝜓 to mean that it is always the case that 𝜑 implies 𝜓 , i.e., for every valuation 𝜎 we have that 𝜎 |= 𝜑 implies 𝜎 |= 𝜓 . Now we describe the semantics of a LinTS. The Semantics of LinTS's. Informally, a LinTS starts at its initial location ℓ * with an arbitrary initial valuation 𝜎 * such that 𝜎 * |= 𝜃 , constituting an initial configuration (ℓ 0 , 𝜎 0 ) = (ℓ * , 𝜎 * ); then at each step 𝑛 (𝑛 ≥ 0), given the current configuration (ℓ 𝑛 , 𝜎 𝑛 ), the LinTS determines the next configuration (ℓ 𝑛+1 , 𝜎 𝑛+1 ) by first selecting a transition 𝜏 = ⟨ℓ, ℓ ′ , 𝜌⟩ such that ℓ = ℓ 𝑛 and then choosing (ℓ 𝑛+1 , 𝜎 𝑛+1 ) to be any configuration that satisfies ℓ 𝑛+1 = ℓ ′ and 𝜎 𝑛 , 𝜎 𝑛+1 |= 𝜌. Formally, the semantics of an LinTS is specified by the notion of paths. A path 𝜋 is a finite sequence of configurations (ℓ 0 , 𝜎 0 ) . . . (ℓ 𝑛 , 𝜎 𝑛 ) such that • (Initialization) ℓ 0 = ℓ * and 𝜎 0 |= 𝜃 , and Intuitively, a path starts with some legitimate initial configuration (as specified by Initialization) and evolves by repeatedly applying the transitions to the current configuration (as described in Consecution). Thus, any path 𝜋 = (ℓ 0 , 𝜎 0 ) . . . (ℓ 𝑛 , 𝜎 𝑛 ) corresponds to a possible evolution of the underlying LinTS.

𝑋 = {𝑥, 𝑦, 𝑡 }, 𝐿 = {ℓ 0 , ℓ 1 }, T = {𝜏 1 , 𝜏 2 }, 𝜏 1 : ⟨ℓ 0 , ℓ 1 , 𝜌 1 ⟩, 𝜏 2 : ⟨ℓ 1 , ℓ 0 , 𝜌 2 ⟩, 𝜃 : 𝑥 = 0 ∧ 𝑦 = 0 ∧ 𝑡 = 0, 𝜌 1 :       𝑡 ′ -𝑡 ≤ 𝑥 ′ -𝑥 ≤ 2(𝑡 ′ -𝑡)∧ 𝑡 ′ -𝑡 ≤ 𝑦 ′ -𝑦 ≤ 2(𝑡 ′ -𝑡)∧ 1 ≤ 𝑡 ′ -𝑡 ≤ 2       , 𝜌 2 :       𝑡 ′ -𝑡 ≤ 𝑥 ′ -𝑥 ≤ 2(𝑡 ′ -𝑡)∧ -(𝑡 ′ -𝑡) ≤ 𝑦 ′ -𝑦 ≤ -2(𝑡 ′ -𝑡)∧ 1 ≤ 𝑡 ′ -𝑡 ≤ 2      
Example 1. Consider a scenario of a vagrant robot from [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF]. The control of the robot works in two alternating modes modelled as two locations ℓ 0 , ℓ 1 . Each mode takes a time between 1 and 2 seconds to complete its task. In mode ℓ 0 , the robot moves in the positive direction of both 𝑥 and 𝑦, and in mode ℓ 1 , it moves in the positive direction of 𝑥 and the negative direction of 𝑦. Figure 1 shows the LinTS of a vagrant robot that consists of three variables 𝑥, 𝑦, 𝑡 and the two locations ℓ 0 , ℓ 1 . The variable pair (𝑥, 𝑦) corresponds to the current position of the robot on the two-dimensional plane, while the variable 𝑡 records the amount of the elapsed time. From the LinTS, we have the following:

• At the start, the robot is at its initial position (0, 0) (i.e., 𝑥 = 0 and 𝑦 = 0) at time 𝑡 = 0 in the initial mode ℓ 0 . • At the location ℓ 0 , the robot takes the transition 𝜏 1 . During the transition, the robot first moves for a time period Δ𝑡 = 𝑡 ′ -𝑡 between 1 to 2 seconds (as indicated by the linear assertion 1 ≤ 𝑡 ′ -𝑡 ≤ 2 in 𝜌 1 ), for which the movement is in the positive direction of both the 𝑥-and 𝑦-axis, each with a nondeterministic distance that falls in the interval [Δ𝑡, 2Δ𝑡] (as specified by the linear assertions Δ𝑡 ≤ 𝑥 ′ -𝑥 ≤ 2Δ𝑡 and Δ𝑡 ≤ 𝑦 ′ -𝑦 ≤ 2Δ𝑡 in 𝜌 1 ); then the robot changes its mode to ℓ 1 . • At the location ℓ 1 , the robot takes the transition 𝜏 2 and changes its mode to ℓ 0 . The movement of the robot is specified by 𝜌 2 and similar to the situation at ℓ 0 . The only difference is that the robot now moves in the positive direction along the 𝑥-axis and in the negative direction along the 𝑦-axis. • The robot switches between ℓ 0 and ℓ 1 by alternately taking the transitions 𝜏 1 and 𝜏 2 .

A path under this LinTS is (ℓ 0 , (𝑥, 𝑦, 𝑡) = (0, 0, 0)), (ℓ 1 , (𝑥, 𝑦, 𝑡) = (2, 2, 1)), (ℓ 0 , (𝑥, 𝑦, 𝑡) = (3, 1, 2)). □

In this work, we consider algorithms for linear invariant generation that work on LinTS's (as the abstract model). Informally, an invariant at a location is a logical formula that is always satisfied by the values of the variables whenever the location is entered by some path. An invariant is linear if it is a linear assertion. The formal definition is as follows.

(Linear) Invariants. An invariant at a location ℓ of a LinTS is a logical formula 𝜑 such that for every path under the LinTS 𝜋 = (ℓ 0 , 𝜎 0 ) . . . (ℓ 𝑛 , 𝜎 𝑛 ) and 0 ≤ 𝑘 ≤ 𝑛, it holds that ℓ 𝑘 = ℓ implies 𝜎 𝑘 |= 𝜑. Furthermore, an invariant 𝜑 is linear if 𝜑 is a linear assertion over the variable set 𝑋 .

To automatically generate invariants, one often investigates a strengthened notion called inductive invariants. Since we only consider linear invariants, we directly present the definition of inductive linear invariants, and in the form of inductive linear assertion maps.

(Inductive) Linear Assertion Maps. A linear assertion map over an LinTS is a function 𝜂 that maps every location ℓ to a linear assertion 𝜂 (ℓ) over the variables 𝑋 . A linear assertion map 𝜂 is inductive if the following conditions hold:

• (Initialization) 𝜃 |= 𝜂 (ℓ * ); • (Consecution) For every transition 𝜏 = ⟨ℓ, ℓ ′ , 𝜌⟩, we have that 𝜂 (ℓ) ∧ 𝜌 |= 𝜂 (ℓ ′ ) ′ , where 𝜂 (ℓ ′ ) ′ is the linear assertion obtained by replacing every variable 𝑥 ∈ 𝑋 in 𝜂 (ℓ ′ ) with its next-value counterpart 𝑥 ′ ∈ 𝑋 ′ .

Informally, a linear assertion map is inductive if it is (i) implied by the initial condition given by 𝜃 at the initial location ℓ * (i.e., Initialization) and (ii) preserved under the application of every transition (i.e., Consecution). By a straightfoward induction on the length of a path under a LinTS, one could verify that the linear assertion at every location in an inductive linear assertion map is guaranteed to be a linear invariant. In the rest of the work, we focus on the automated synthesis of inductive linear assertion maps.

AN OVERVIEW OF OUR APPROACH

In this section, we first review the original approaches in [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that generate linear invariants through Farkas' Lemma and point out the weakness of each approach, and then sketch our key improvements to these approaches.

3.1 The Original Approaches [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] In [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], a sound and complete constraint-solving framework for linear invariant generation is proposed via Farkas' Lemma [START_REF] Farkas | A Fourier-féle mechanikai elv alkalmazásai (Hungarian)[END_REF]. The use of Farkas' Lemma transforms the inductive condition for linear invariants equivalently into a system of quadratic constraints, by solving which one could obtain concrete linear invariants. The key merit of the use of Farkas' Lemma is that it simplifies the constraints from the inductive condition to quadratic constraints. In [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF], the quadratic constraints obtained after applying Farkas' Lemma were solved exactly through quantifier elimination. While in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], these constraints were solved instead through polyhedra manipulation by applying several reasonable heuristics to avoid the high runtime complexity caused by quantifier elimination.

To review these two approaches, we first recall Farkas' Lemma. Farkas' Lemma is a fundamental theorem that characterizes basic relationships between linear inequalities. Below we present the version of Farkas' Lemma that charactertizes the entailment from a linear assertion to a linear inequality. We follow the presentation form of Farkas' Lemma in [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF]. Theorem 3.1 (Farkas' Lemma). Consider a linear assertion 𝜑 over a set 𝑉 = {𝑥 1 , . . . , 𝑥 𝑛 } of real-valued variables in the form of a conjunction of the following linear inequalities: 𝜑 :

𝑎 11 • 𝑥 1 + • • • + 𝑎 1𝑛 • 𝑥 𝑛 + 𝑏 1 ≥ 0 . . . . . . . . . 𝑎 𝑚1 • 𝑥 1 + • • • + 𝑎 𝑚𝑛 • 𝑥 𝑛 + 𝑏 𝑚 ≥ 0
When 𝜑 is satisfiable (i.e., there is a valuation over 𝑉 that satisfies 𝜑), it implies a linear inequality 𝜓 𝜓 :

𝑐 1 • 𝑥 1 + • • • + 𝑐 𝑛 • 𝑥 𝑛 + 𝑑 ≥ 0 (i.
𝜆 0 1 ≥ 0 𝜆 1 𝑎 11 • 𝑥 1 + • • • + 𝑎 1𝑛 • 𝑥 𝑛 + 𝑏 1 1 0 . . . . . . . . . . . . 𝜆 𝑚 𝑎 𝑚1 • 𝑥 1 + • • • + 𝑎 𝑚𝑛 • 𝑥 𝑛 + 𝑏 𝑚 𝑚 0 𝑐 1 • 𝑥 1 + • • • + 𝑐 𝑛 • 𝑥 𝑛 + 𝑑 ≥ 0 -1 ≥ 0          𝜑 ← 𝜓 ← false
One direction of Farkas' Lemma is straightforward, as one easily sees that if we have a nonnegative linear combination of the inequalities in 𝜑 that can derive 𝜓 , then it is guaranteed that 𝜓 holds whenever 𝜑 is true. Farkas' Lemma further establishes that the other direction is also valid. Remark 1. In the statement of Farkas' Lemma above, if we change a linear inequality 𝑎 𝑗1 𝑥 1 + • • • + 𝑎 𝑗𝑛 𝑥 𝑛 + 𝑏 𝑗 ≥ 0 in 𝜑 to equality (i.e., 𝑎 𝑗1 𝑥 1 + • • • + 𝑎 𝑗𝑛 𝑥 𝑛 + 𝑏 𝑗 = 0), then the theorem holds with the relaxation that we do not require 𝜆 𝑗 ≥ 0. This could be easily observed by first replacing the equality equivalent with both 𝑎 𝑗1 𝑥 1 + • • • + 𝑎 𝑗𝑛 𝑥 𝑛 + 𝑏 𝑗 ≥ 0 and 𝑎 𝑗1 𝑥 1 + • • • + 𝑎 𝑗𝑛 𝑥 𝑛 + 𝑏 𝑗 ≤ 0, and then applying Farkas' Lemma. By similar arguments, the theorem statement holds upon changing multiple linear inequalities into equalities with the relaxation of non-negativity for their corresponding 𝜆 𝑗 's.

The application of Farkas' Lemma can be visualized by the tabular form in Table 1, where 1 , . . . , 𝑚 ∈ {=, ≥} and we mulitply 𝜆 0 , 𝜆 1 , . . . , 𝜆 𝑚 with their inequalities in 𝜑 and sum up them together to get 𝜓 . For 1 ≤ 𝑗 ≤ 𝑚, if 𝑗 is ≥, we require 𝜆 𝑗 ≥ 0, otherwise (i.e., 𝑗 is =) we do not impose restriction on 𝜆 𝑗 . We then recall several concepts from polyhedra theory. Polyhedra and polyhedral cones. A subset 𝑃 of R 𝑛 is a polyhedron if 𝑃 = {x ∈ R 𝑛 | A • x ≤ b} for some real matrix 𝐴 ∈ R 𝑚×𝑛 and real vector b ∈ R 𝑚 , where x is treated as a column vector and the comparison A • x ≤ b is defined in the coordinate-wise fashion. A polyhedron 𝑃 is a polyhedral cone if 𝑃 = {x ∈ R 𝑛 | A • x ≤ 0} for some real matrix 𝐴 ∈ R 𝑚×𝑛 , where 0 is the 𝑚-dimensional zero column vector. It is well-known from Farkas-Minkowski-Weyl Theorem [START_REF] Schrijver | Theory of linear and integer programming[END_REF]Corollary 7.1a] that any polyhedral cone 𝑃 can be represented as 𝑃 = { 𝑘 𝑖=1 𝜆 𝑖 • g 𝑖 | 𝜆 𝑖 ≥ 0 for all 1 ≤ 𝑖 ≤ 𝑘 } for some real vectors g 1 , . . . , g 𝑘 , where such real vectors g 𝑖 's are called a collection of generators for the polyhedral cone 𝑃.

Polyhedron projection. For a polyhedron 𝑃

= {(x T , u T ) T ∈ R 𝑝+𝑞 | A • x + B • u ≤ c} where A ∈ R 𝑚×𝑝 , B ∈ R 𝑚×𝑞
are real matrices and 𝑐 ∈ R 𝑚 is a real vector, the projection of 𝑃 onto the dimensions x (i.e., the first 𝑝 dimensions) is defined as the polyhedron 𝑃

[x] := {x ∈ R 𝑝 | ∃u ∈ R 𝑞 .(x T , u T ) T ∈ 𝑃 }.
It is guaranteed by e.g. Fourier-Motzkin Elimination [64, Chapter 12.2] that 𝑃 [x] will always be a polyhedron. Now we review the approaches in [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. The invariant-generation workflow. Based on Farkas' Lemma, the approaches in [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] generate linear invariants over a LinTS by the following steps (Steps A1 -A3). Below we fix an input LinTS with variables 𝑋 = {𝑥 1 , . . . , 𝑥 𝑛 }.

Step A1 In the first step, both the approaches establish a template for an inductive linear assertion maps. A template 𝜂 involves a linear inequality

𝜂 (ℓ) = 𝑐 ℓ,1 • 𝑥 1 + • • • +𝑐 ℓ,𝑛 • 𝑥 𝑛 +𝑑 ≥ 0 at
each location ℓ of the LinTS, such that the coefficients 𝑐 ℓ,1 , . . . , 𝑐 ℓ,𝑛 , 𝑑 are unknown and to be resolved.

Step A2 In the second step, both the approaches establish constraints from the initialization and the consecution conditions for invariants. The initialization condition specifies that the linear inequality 𝜂 (ℓ * ) at the initial location ℓ * should be implied by the initial condition 𝜃 , i.e., 𝜃 |= 𝜂 (ℓ * ). The consecution condition specifies that every transition preserves the linear assertion map 𝜂, i.e., for every transition ⟨ℓ, ℓ ′ , 𝜌⟩ we have that 𝜂 (ℓ) ∧ 𝜌 |= 𝜂 (ℓ ′ ) ′ .

Step A3 In the third step, both the approaches apply Farkas' Lemma to the constraints collected from the initialization condition 𝜃 |= 𝜂 (ℓ * ) and the consecution conditions 𝜂 (ℓ) ∧ 𝜌 |= 𝜂 (ℓ ′ ) ′ for each transition ⟨ℓ, ℓ ′ , 𝜌⟩. For intialization, we apply the tabular form (Table 1) to obtain

𝜆 0 1 ≥ 0 𝜆 1 𝑎 11 𝑥 1 + • • • + 𝑎 1𝑛 𝑥 𝑛 + 𝑏 1 1 0 . . . . . . . . . . . . 𝜆 𝑚 𝑎 𝑚1 𝑥 1 + • • • + 𝑎 𝑚𝑛 𝑥 𝑛 + 𝑏 𝑚 𝑚 0 𝑐 ℓ * ,1 𝑥 1 + • • • + 𝑐 ℓ * ,𝑛 𝑥 𝑛 + 𝑑 ℓ * ≥ 0 -1 ≥ 0          𝜃 ← 𝜂 (ℓ * ) ← false
which results in a linear assertion over the unknown coefficients 𝑐 ℓ * ,1 , . . . , 𝑐 ℓ * ,𝑛 , 𝑑 and the fresh variables 𝜆 0 , 𝜆 1 , . . . , 𝜆 𝑚 . Similarly, the tabular form for the consecution condition of a transition ⟨ℓ, ℓ ′ , 𝜌⟩ gives

𝜇 𝑐 ℓ,1 𝑥 1 + • • • + 𝑐 ℓ,𝑛 𝑥 𝑛 + 𝑑 ℓ ≥ 0 𝜆 0 1 ≥ 0 𝜆 1 𝑎 11 𝑥 1 + • • • + 𝑎 1𝑛 𝑥 𝑛 + 𝑎 ′ 11 𝑥 ′ 1 + • • • + 𝑎 ′ 1𝑛 𝑥 ′ 𝑛 + 𝑏 1 1 0 . . . . . . . . . . . . . . . . . . 𝜆 𝑚 𝑎 𝑚1 𝑥 1 + • • • +𝑎 𝑚𝑛 𝑥 𝑛 +𝑎 ′ 𝑚1 𝑥 ′ 1 + • • • +𝑎 ′ 𝑚𝑛 𝑥 ′ 𝑛 +𝑏 𝑚 𝑚 0 𝑐 ℓ ′ ,1 𝑥 ′ 1 + • • • +𝑐 ℓ ′ ,𝑛 𝑥 ′ 𝑛 + 𝑑 ℓ ′ ≥ 0 -1 ≥ 0 ← 𝜂 (ℓ)          𝜌 ←𝜂 (ℓ ′ ) ′ ← false
where in addition to 𝜆 𝑗 , 𝑐 ℓ,𝑗 , 𝑑 ℓ , 𝑐 ℓ ′ ,𝑗 , 𝑑 ℓ ′ we have a fresh variable 𝜇 as the non-negative muliplier for 𝜂 (ℓ). Note that for the consecution condition, the constraint obtained is no longer linear since the fresh variable 𝜇 is multiplied to 𝜂 (ℓ) in the tabular above.

Step A4 In the last step, the (non-linear) constraints collected in the previous step are solved to obtain the concrete values for the unknown coefficients in 𝜂, so that a concrete inductive linear assertion map would be obtained. In [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF], these constraints were solved through quantifier elimination, while in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] the constraints were solved through (i) several reasonable heuristics to guess possible values for the key parameter 𝜇 so as to remove the non-linearity and (ii) the generator computation of a polyhedral cone originally represented in its linear inequalities. A major heuristics adopted in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] to guess possible values for 𝜇 is based on some practical rules such as factorization and setting 𝜇 manually to 0, 1 (where 0 means an invariant local to the guard of the transition and 1 means an invariant incremental to the previous step).

Below we focus on the workflow of [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], since it is the most related work to our result. Since the approach of [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] resolves the mulipliers 𝜇 for consecution by guessing their concrete values, the guessed values for the consecution condition of every transition in a LinTS result in a disjunction of linear assertions over the unknown coefficients in the template 𝜂, where each guessed value for the parameter 𝜇 corresponds to one linear assertion in the disjunction. By combining conjunctively the constraints obtained from the consecution condition (for every tranisition) and those from the initial condition, we obtain a propositional formula in conjunctive normal form (CNF) where each atomic proposition is a linear assertion over the unknown coefficients of 𝜂. More specifically, each atomic proposition is in the form A • c ≤ 0 where c is the vector of unknown coefficients in the template, hence is a polyhedral cone. Having obtained the CNF formula, the approach expands the formula (through the distributive law between conjunction and disjunction) equivalently into a DNF (disjunctive normal form) formula where each disjunctive clause is a conjunction of certain atomic propositions in the original CNF formula (therefore being still a polyhedral cone), and then computes through the double description method [START_REF] Bagnara | Possibly Not Closed Convex Polyhedra and the Parma Polyhedra Library[END_REF] and Parma Polyhedra library (PPL) [START_REF] Bagnara | Possibly Not Closed Convex Polyhedra and the Parma Polyhedra Library[END_REF][START_REF] Ppl | [END_REF] the generators of each disjunctive clause (in the form of a polyhedral cone) once a disjunctive clause is expanded. By instantiating the computed generators back to the unknown coefficients in the template, the approach obtains the solved concrete linear invariants, for which one generator corresponds to one inductive linear asertion map.

A key improvement to reduce the runtime arising from the expansion from CNF to DNF in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] is the subsumption testing that checks whether the current conjunction of atomic propositions during the expansion is already subsumed by (i.e., implies) the invariants obtained from the previouslyexpanded disjunctive clauses; once the subsumption is successful, one knows that the current conjunction will not produce more meaningful invariants, hence the computation for the current conjunction halts and the expansion proceeds by switching to other branches. Below we present an example to illustrate the workflow of [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF].

Example 2. Consider the LinTS in Example 1. The approach in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] first establishes a template 𝜂 by setting 𝜂 (ℓ

𝑖 ) := 𝑐 ℓ 𝑖 ,1 𝑥 + 𝑐 ℓ 𝑖 ,2 𝑦 + 𝑐 ℓ 𝑖 ,3 𝑡 + 𝑑 ℓ 𝑖 ≥ 0 for 𝑖 ∈ {0, 1}.
Then the approach encodes initialization and consecution by the tabular form in Table 1. For the initialization, we have

𝜆 0 1 ≥ 0 𝜆 1 𝑥 = 0 𝜆 2 𝑦 = 0 𝜆 3 𝑡 = 0 𝑐 ℓ 0 ,1 𝑥 + 𝑐 ℓ 0 ,2 𝑦 + 𝑐 ℓ 0 ,3 𝑡 + 𝑑 ℓ 0 ≥ 0        𝜃 ← 𝜂 (ℓ 0 )
that results in the constraints 𝑐 ℓ 0 ,1 = 𝜆 1 , 𝑐 ℓ 0 ,2 = 𝜆 2 , 𝑐 ℓ 0 ,3 = 𝜆 3 , 𝑑 ℓ 0 ≥ 0 . After projecting away the fresh variables 𝜆 𝑗 's, we obtain [𝑑 ℓ 0 ≥ 0] for the initialization. For the consecution conditions, we present the tabular form for the transition 𝜏 1 :

𝜇 𝑐 ℓ 0 ,1 𝑥 +𝑐 ℓ 0 ,2 𝑦 +𝑐 ℓ 0 ,3 𝑡 +𝑑 ℓ 0 ≥ 0 𝜆 0 1 ≥ 0 𝜆 1 -𝑥 + 𝑡 + 𝑥 ′ - 𝑡 ′ ≥ 0 𝜆 2 𝑥 -2𝑡 - 𝑥 ′ + 2𝑡 ′ ≥ 0 𝜆 3 -𝑦 -2𝑡 + 𝑦 ′ + 2𝑡 ′ ≥ 0 𝜆 4 𝑦 + 𝑡 - 𝑦 ′ - 𝑡 ′ ≥ 0 𝜆 5 - 𝑡 𝑡 ′ -1 ≥ 0 𝜆 6 𝑡 - 𝑡 ′ + 2 ≥ 0 𝑐 ℓ 1 ,1 𝑥 ′ +𝑐 ℓ 1 ,2 𝑦 ′ +𝑐 ℓ 1 ,3 𝑡 ′ +𝑑 ℓ 1 ≥ 0 ← 𝜂 (ℓ 0 )                    𝜌 1 ←𝜂 (ℓ 1 ) ′
where the fresh variables 𝜆 𝑗 's are local to the tabular above and does not overlap with the tabulars for other transitions. The 𝜆 𝑗 's are again eliminated by polyhedron projection, while the fresh variable 𝜇 is eliminated by heuristically guessing its values. The transition 𝜏 2 is treated in the similar way.

After guessing the values for 𝜇 with several heuristics, the resultant constraint is a CNF formula where each atomic proposition is a polyhedral cone over the unknown coefficients. The approach further expands the CNF equivalently into a DNF formula (which we omit here due to its vast amount of technicality). For this example, one disjunctive clause from the DNF formula is as follows (where we abbreviate 𝑐 ℓ 𝑖 ,𝑗 as 𝑐 𝑖 𝑗 , 𝑑 ℓ 𝑖 as 𝑑 𝑖 ):

            𝑐 01 =𝑐 11 , 𝑐 02 =𝑐 12 , 𝑐 03 =𝑐 13 , 𝑑 0 ≥ 0, 2𝑐 01 -𝑑 0 +𝑐 12 +𝑐 13 +𝑑 1 ≥ 0, 2𝑐 01 -𝑑 0 +2𝑐 12 +𝑐 13 +𝑑 1 ≥ 0, 4𝑐 01 -𝑑 0 +4𝑐 12 +2𝑐 13 +𝑑 1 ≥ 0, 𝑐 01 -𝑑 0 +𝑐 12 +𝑐 13 +𝑑 1 ≥ 0, 𝑐 01 -𝑑 0 +2𝑐 12 +𝑐 13 +𝑑 1 ≥ 0,2𝑐 01 -𝑑 0 +4𝑐 12 +2𝑐 13 +𝑑 1 ≥ 0, 2𝑐 01 +𝑑 0 -4𝑐 12 +2𝑐 13 -𝑑 1 ≥ 0,4𝑐 01 +𝑑 0 -4𝑐 12 +2𝑐 13 -𝑑 1 ≥ 0, 𝑐 01 +𝑑 0 -2𝑐 12 +𝑐 13 -𝑑 1 ≥ 0, 2𝑐 01 +𝑑 0 -2𝑐 12 +𝑐 13 -𝑑 1 ≥ 0, 𝑐 01 +𝑑 0 -𝑐 12 +𝑐 13 -𝑑 1 ≥ 0, 2𝑐 01 +𝑑 0 -𝑐 12 +𝑐 13 -𝑑 1 ≥ 0             (1)
By computing the generators of the polyhedral cone above, we obtain the following generators and their corresponding invariants in Table 2, where in the left part each row specifies a generator (over the unknown coefficients 𝑐 𝑖 𝑗 , 𝑑 𝑖 's) and in the right part we instantiate the generator to the unknown coefficients in the template 𝜂 to obtain the invariants at locations ℓ 0 and ℓ 1 .

Table 2. Generators (left) and Their Invariants (right) for ( 1)

𝑐 01 𝑐 02 𝑐 03 𝑑 0 𝑐 11 𝑐 12 𝑐 13 𝑑 1 𝜂 (ℓ 0 ) 𝜂 (ℓ 1 ) 0 0 0 1 0 0 0 1 1 ≥ 0 1 ≥ 0 0 0 0 0 0 0 0 0 0 ≥ 0 0 ≥ 0 -1 0 2 0 -1 0 2 0 -𝑥 + 2𝑡 ≥ 0 -𝑥 + 2𝑡 ≥ 0 0 -1 1 0 0 -1 1 2 -𝑦 + 𝑡 ≥ 0 -𝑦 + 𝑡 + 2 ≥ 0 0 -1 2 0 0 -1 2 0 -𝑦 + 2𝑡 ≥ 0 -𝑦 + 2𝑡 ≥ 0 0 0 1 0 0 0 1 1 𝑡 ≥ 0 𝑡 + 1 ≥ 0 0 1 2 0 0 1 2 0 𝑦 + 2𝑡 ≥ 0 𝑦 + 2𝑡 ≥ 0 0 1 1 0 0 1 1 -2 𝑦 + 𝑡 ≥ 0 𝑦 + 𝑡 -2 ≥ 0 1 0 -1 0 1 0 -1 0 𝑥 -𝑡 ≥ 0 𝑥 -𝑡 ≥ 0 0 0 1 0 0 0 1 -1 𝑡 ≥ 0 𝑡 -1 ≥ 0
The obtained invariants are further minimized, so that the final invariants obtained at 𝜂 (ℓ 0 ) are

-𝑥 + 2𝑡 ≥ 0 -𝑦 + 𝑡 ≥ 0 𝑥 -𝑡 ≥ 0 𝑦 + 𝑡 ≥ 0
and for 𝜂 (ℓ 1 ) the invariants are

-𝑥 + 2𝑡 ≥ 0 -𝑦 + 𝑡 + 2 ≥ 0 -𝑦 + 2𝑡 ≥ 0 𝑡 -1 ≥ 0 𝑦 + 𝑡 -2 ≥ 0 𝑥 -𝑡 ≥ 0 .
It happens that the invariants from the whole DNF formula coincide with those computed from (1). □

Weakness of [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. Although the approaches [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] provides the first and only invariant generation approaches via Farkas' Lemma, they have the following weakness. The approach [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF] generates the invariants by quantifier elimination, thus is impractical [START_REF] Zhao | A "Hybrid" Approach for Synthesizing Optimal Controllers of Hybrid Systems: A Case Study of the Oil Pump Industrial Example[END_REF]. The approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] avoids the high runtime complexity from quantifier elimination by several reasonable heuristics in the application of Farkas' Lemma, solves the invariants by generator computation, and performs subsumption testing to reduce the combinatorial explosion from the CNF-to-DNF expansion; however, the CNF-to-DNF expansion with subsumption testing still leads to a large amount of combinatorial explosion, and the generator computation also causes a considerable amount of combinatorial explosion.

In order to further reduce combinatorial explosion in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], we complement the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] with (i) a location-by-location idea that generates the invariants one program location at a time and (ii) a segmented subsumption testing that is independent of the original subsumption testing in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. The following subsections illustrate the main idea of our improvements.

Location-by-Location Linear Invariant Generation

A drawback of the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] is that in its last step (Step A4 in the previous subsection), the invariants are obtained as the generators of a whole polyhedral cone that involves the unknown coefficients at all locations. This leads to two obstacles that may largely increase the runtime. The first is that in the subsumption testing, there is only a small chance that a whole polyhedral cone with all the unknown coefficients be subsumed by the invariants already generated. The second is that the generator computation of a polyhedral cone is an expensive operation since it may cause exponential blowup in the number of its variables, thus to compute the generators of a whole polyhedral cone that involves all the unknown coefficients may induce a large amount of runtime.

To address the two obstacles, we propose the novel idea of generating the invariants one location at a time (i.e., location-by-location). With this idea, the chance of successful subsumption in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] is increased since now the subsumption testing involves only one location. Moreover, this idea directly leads to a parallel processing that breaks down the whole invariant generation task at all locations into multiple processors for which each processor only handles a small number of program locations. Note that although the idea seems simple, it can directly lead to improvement on the runtime, and enables further technical improvements. Below we present an example to illustrate our idea.

Example 3. Consider Example 2 again. We run Steps A1-A3 of the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] and the part of Step A4 that computes a CNF from the previous steps. Recall that originally in Step A4, the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] expands the CNF into its equivalent DNF, and computes the generators for each atomic proposition (polyhedral cone) of the DNF. Our location-by-location idea distinguish itself by focusing on one target location in the whole invariant generation process and generate only the invariants at the target location accordingly, and a simple way to implement this idea is to project the polyhedral cones in the DNF onto the unknown coefficients related to the target location. For example, suppose that we focus on the target location ℓ 0 . Then a simple way to implement our idea is to project the polyhedral cones in the DNF onto the unknown coefficients related to ℓ 0 . For instance, for the polyhedral cone in 1, we project the polyhedral cone onto the dimensions specified by the unknown coefficients 𝑐 0𝑗 's (1 ≤ 𝑗 ≤ 3) and 𝑑 0 that are related to the location ℓ 0 , to obtain the following polyhedral cone

𝑑 0 ≥ 0 𝑐 01 -𝑐 02 + 𝑐 03 ≥ 0 𝑐 01 + 𝑐 02 + 𝑐 03 ≥ 0 2𝑐 01 -𝑐 02 + 𝑐 03 ≥ 0 2𝑐 01 + 𝑐 02 + 𝑐 03 ≥ 0 .
The generators of this projected polyhedral cone gives the following linear invariants at the location ℓ 0 :

𝑐 01 𝑐 02 𝑐 03 𝑑 0 𝜂 (ℓ 0 ) 0 0 0 1 1 ≥ 0 0 0 0 0 0 ≥ 0 1 0 -1 0 𝑥 -𝑡 ≥ 0 0 1 1 0 𝑦 + 𝑡 ≥ 0 -1 0 2 0 -𝑥 + 2𝑡 ≥ 0 0 -1 1 0 -𝑦 + 𝑡 ≥ 0
After minimization, the invariants added for the location ℓ 0 from the polyhedral cone (1) include

-𝑥 + 2𝑡 ≥ 0 -𝑦 + 𝑡 ≥ 0 𝑥 -𝑡 ≥ 0 𝑦 + 𝑡 ≥ 0 .
which happen to include all the linear inductive invariants at ℓ 0 . Note that under our idea, we only generate the invariants at ℓ 0 in one invariant-generation process, and keep the invariants at other locations (i.e., ℓ 1 ) to be true. By another invariant-generation process for the location ℓ 1 , we get the invariants at the location ℓ 1 as follows:

-𝑥 + 2𝑡 ≥ 0-𝑦 + 𝑡 + 2 ≥ 0-𝑦 + 2𝑡 ≥ 0 𝑡 -1 ≥ 0 𝑦 + 𝑡 -2 ≥ 0 𝑥 -𝑡 ≥ 0 .
The invariants obtained coincide with the approach in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. □

Segmented Subsumption Testing

We also propose a novel subsumption testing on the CNF-to-DNF expansion in Step A4 of Section 3.1.

The main idea is that the subsumption testing divides the CNF formula into several segments and perform local subsumption testing in each segment. Our subsumption testing is independent of the original one in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] since we focus on local subsumption in each segment, while the original one [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] focuses however on global subsumption testing on the whole CNF formula. Thus, our subsumption testing is able to find local subsumptions in advance so that the original CNF formula can be simplified before it is used to generate invariants.

Below we present a scenario to illustrate this idea.

Example 4. Consider the scenario that the CNF formula to be expanded in

Step A4 is Φ = 𝐶 1 ∧ 𝐶 2 ∧ 𝐶 3 ∧ 𝐶 4
where each 𝐶 𝑖 is a disjunction of polyhedral cones. The idea of segmented subsumption divides the CNF into two segments 𝐶 1 ∧ 𝐶 2 and 𝐶 3 ∧ 𝐶 4 , and perform local subsumption testing in each segment. For the segment 𝐶 1 ∧ 𝐶 2 (as a sub-CNF formula), we expand it into a (sub-)DNF and remove redundant disjunctive clauses by subsumption testing, to obtain a DNF formula Θ 1 . We also perform the same operation to 𝐶 3 ∧ 𝐶 4 and obtain a DNF Θ 2 . Finally, we apply our improvements in the previous subsection to the CNF Θ 1 ∧ Θ 2 (which is equivalent to Φ) to generate the invariants. □

LOCATION-BY-LOCATION LINEAR INVARIANT GENERATION

In this section, we formally demonstrate our idea of generating the linear invariants at only one target program location in each invariant generation process, thus generating the linear invariants location-by-location. The idea could be easily implemented by projecting the generated invariants onto the (related coefficents at) the target location in Step A4 of Section 3.1, and leads directly to a parallel processing that assigns the invariant generation task at every program location separately onto different processors. Based on the idea, in the following we present two technical improvements that further reinforce the idea by reducing more combinatorial explosion in the original approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. The first improvement is a reordering of the expansion from the CNF into its equivalent DNF at Step A4 (Section 3.1) so that the transitions associated with the target program location are expanded first to detect possible subsumptions earlier. The second improvement is a technique to compute the generators (as the generated invariants, see Step A4 of Section 3.1) at the target location via the projection onto the location. Below we fix a LinTS and a target location ℓ ★ on which invariant generation is considered.

Reordered CNF-to-DNF Expansion

Recall that in Step A4, we obtain a CNF formula and expand the CNF into its equivalent DNF from which we generate invariants. Since the transformation from CNF into its DNF causes combinatorial explosion, the original approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] adopts a subsumption testing that detects redundant disjunctive clauses in the DNF. Here, since we focus only on the target location ℓ ★ , we can reorder the expansion for the CNF formula by expanding the transitions mostly related to the target location first so that subsumptions could be detected much earlier than the original approach. The technical details is as follows. Reordering of the CNF-to-DNF Expansion. We reorder the expansion from the CNF to its DNF so that the transitions of the LinTS associated with ℓ ★ are expanded first, and the other transitions are expanded last. The reason to expand the transitions with the target location ℓ ★ first is as follows: since we focus on the invariant generation at ℓ ★ , successful subsumptions would be detected earlier if we check transitions that involve the target location ℓ ★ first; in contrast, if one chooses an arbitary order for the expansion, then from our idea the expansion at locations other than ℓ ★ cannot lead to successful subsumptions, as the subsumption testing is operated by checking the polyhedral inclusion that only involves the unknown coefficients at ℓ ★ (i.e., not involving the unknown coefficients at other locations). We further refine the epansion order by having that the intra-transitions in the LinTS that involves only the location ℓ ★ (i.e., the transitions from ℓ ★ to itself) are expanded first, those that involves ℓ ★ and other locations (i.e., the inter-transitions from ℓ ★ to other locations or from other locations to ℓ ★ ) are expanded second, those that does not involve location ℓ ★ are expanded last; the heuristics behind the refined reordering is that they are ordered by the descending order of relevance to the target location ℓ ★ .

The pseudo-code for the invariant-generation algorithm that adopts the location-by-location idea and integrates the reordered expansion is given in Algorithm 1. In the pseudo-code, the algorithm first reorders the CNF formula 𝑖 𝐶 𝑖 obtained from Step A4 in Section 3.1 as stated in the previous paragraph (line 1). The reordered CNF is 𝑖 𝐶 ′ 𝑖 . Then the algorithm initializes the values for inv(ℓ ★ ), 𝑖 and all 𝑛 𝑖 's (lines 2-4): the variable inv(ℓ ★ ) is a set variable used to collect linear invariants generated in the algorithm and initialized to be ∅; since we implement expansion from the reordered CNF to its DNF through backtracking, we use the variable 𝑖 to represent that the algorithm is currently traversing the atomic propositions in 𝐶 ′ 𝑖 , and the variable 𝑛 𝑖 to represent that the current atomic proposition being traversed in 𝐶 ′ 𝑖 is the 𝑛 𝑖 -th atomic proposition. Next, the algorithm follows the backtracking process that iteratively selects one atomic proposition 𝑑 𝑖 from each 𝐶 ′ 𝑖 (line 7) following the increasing order of 𝑖, aggressively checks whether the current exploration has already been subsumed by the current invariant inv(ℓ ★ ) (line 9), and updates the current invariant set by adding the invariants that can be derived from the conjunction 𝑚 𝑠=1 𝑑 𝑠 (as a polyhedral cone) and only involve the unknown coefficients at the target location ℓ ★ , which is done through the Gen procedure if all the preceding subsumption testings have been unsuccessful. Note that in the subsumption testing, the condition 𝑖 𝑠=1 𝑑 𝑠 ⊈ inv(ℓ ★ ) means that the polyhedral cone defined by the linear assertion 𝑖 𝑠=1 𝑑 𝑠 is not a subset of the counterpart generated by the linear invariants (as generators extended to full dimension) in inv(ℓ ★ ). The return value is the final set value of inv(ℓ ★ ) after the whole backtracking process, which collects all the linear invariants that are generators of the polyhedral cones having passed the subsumption testing.

A missing part in the pseudo-code is the details for the procedure Gen, whose counterpart in the original work [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] is the generator computation of the polyhedral cone that involves the unknown coefficients at all program locations. In the next part, we instantiate the procedure by our second improvement related to location-by-location linear invariant generation that proposes various procedures to compute the generators that only involve the dimensions at the target location.

Projection-Based Generator Computation

In the original work [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], the final invariants are generated by computing the generators of the polyhedral cones (i.e., poly at line 11 in Algorithm 1) obtained from the transformed DNF. As stated previously, the polyhedral cones involve the unknown coefficients at all program locations, thus the generation computation over these polyhedral cones may cause combinatorial explosion. From our idea of generating the invariants only at the target location, one only needs to compute the generators that involve the dimensions at the target location. This allows further technical Algorithm 1 Linear Invariant Generation at Location ℓ ★ Input: ℓ ★ : target location; 𝑚 𝑖=1 𝐶 𝑖 : the CNF formula obtained from Step A4 (Section 3.1) in the approach [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF] where each 𝐶 𝑖 is the 𝑖-th conjunctive clause and is a disjunction of polyhedra over the unknown coefficients in the template; Reordering( 𝑖 𝐶 𝑖 ) : a procedure that reorders all 𝐶 𝑖 's into 𝑖 𝐶 ′ 𝑖 w.r.t whether each 𝐶 𝑖 corresponds to an intra-or an inter-transition as stated previously; we write

𝐶 ′ 𝑖 = 𝑁 𝑖 𝑗=1 𝑐 ′ 𝑖,𝑗
where each 𝑐 ′ 𝑖,𝑗 is an atomic proposition in the form of a polyhedron over the unknown coefficients; Gen(ℓ, poly) : a procedure that generates the invariants at a location ℓ given a polyhedron poly over the unknown coefficients Output: a collection inv(ℓ ★ ) of linear invariants at ℓ ★ ; end if 20: end while; 21: return inv(ℓ ★ ); improvements that can reduce the combinatorial explosion from the generator computation. Below we propose several solutions that compute generators only at the dimensions of the target location. These solutions are to be substituted into the procedure Gen in the pseudo-code of Algorithm 1 to complete our invariant-generation procedure for location-by-location linear invariant generation. Generator Computation at Target Location We propose four solutions to compute the invariants at the dimensions of the target location, namely generator projection, Fourier-Motzkin Elimination (FME), Block Elimination and our improvement on Block Elimination. Below we fix a polyhedral cone 𝑃 = {(x T , u T ) T ∈ R 𝑝+𝑞 | Ax + Bu ≤ c} and consider to compute the generators only over the dimensions x that correspond to the dimensions of unknown coefficients at our target location ℓ ★ . Generator Projection. Our first solution is still to compute the generators directly for the polyhedral cones that involves all program locations by standard methods (e.g., the double-description method [START_REF] Fukuda | Double Description Method Revisited[END_REF]), but to project the computed generators onto the unknown coefficients at the target location ℓ ★ . Thus in detail, the first solution computes the generators of the polyhedral cone 𝑃 and project the generators to the dimensions x. Fourier-Motzkin Elimination. Our second solution is to apply the classical method of Fourier-Motzkin Elimination (FME) [START_REF] Schrijver | Theory of linear and integer programming[END_REF]Chapter 12.2] that first projects away the extra dimensions other than the unknown coefficients at the target location ℓ ★ and then computes the generators. In detail, the solution through FME first eliminates the variables u one by one in 𝑃 to get the projected polyhedral cone 𝑃 [x] on the dimensions x, then compute the generators of 𝑃 [x] by standard methods. Block Elimination. We consider applying Block Elimination (see e.g. [START_REF] Talaashrafi | Complexity Results for Fourier-Motzkin Elimination[END_REF]) as our third solution, which is an alternative approach that eliminate unrelated variables through the computation of the generators of the projection cone of a polyhedron. The motivation is to address the issue of generating possibly many redundant inequalities in the application of FME. The details of Block Elimination is as follows. To illustrate this method, we first present a variant form of Farkas' Lemma [START_REF] Schrijver | Theory of linear and integer programming[END_REF]Corollary 7.1e]. Theorem 4.1. Let A be a real matrix and b be a real vector. The polyhedron 𝑃 = {x | Ax ≤ b} is non-empty if and only if there is no vector y satisfying y ≥ 0, y T A = 0 and y T b < 0.

The polyhedral cone {y | y ≥ 0 ∧ y T A = 0} derived from the polyhedron 𝑃 = {x | Ax ≤ b} in the statement of Theorem 4.1 is usually called the projection cone of 𝑃, which we denoted by proj(A) since it is only related to the matrix A (that help defines 𝑃).

Block Elimination utilizes projection cones to reduce the amount of redundant inequalities that may arise in FME. Consider the polyhedron 𝑃 = {(x T , u T ) T ∈ R 𝑝+𝑞 | Ax + Bu ≤ c} and the projection of 𝑃 onto x. The method treats 𝑃 as a polyhedron 𝑄 (x) with the parameter x such that 𝑄 (x) := {u ∈ R 𝑞 | Bu ≤ c -Ax} and computes the generators of the projection cone proj(B). By Theorem 4.1, to ensure that 𝑄 (x) is non-empty, it suffices to guarantee that g T (c -Ax) ≥ 0 for all the generators g of proj(𝑄 (x)). Hence, after the generators g 1 , . . . , g 𝑘 of the projection cone proj(𝑄 (x)) are computed, the method outputs the projected polyhedron 𝑃 [x] as defined by the linear inequalities g T 𝑗 (c -Ax) ≥ 0 for 1 ≤ 𝑗 ≤ 𝑘. The remaining task is again to compute the generators of the projected polyhedral cone 𝑃 [x]. Our improvement on Block Elimination. Recall that our task is to compute the projected generators of a polyhedral cone Ac ≤ 0 onto the unknown coefficients at ℓ ★ (here c is the vector of all unknown coefficients in the template), where the polyhedral cone Ac ≤ 0 is a disjunctive clause from the DNF obtained in Step A4. Following Block Elimination, we still project the polyhedral cone Ac ≤ 0 onto the unknown coefficients at ℓ ★ , but have further improvements as follows.

Below we describe our improvements in detail. Let c ★ be the vector of unknown coefficients at ℓ ★ and c * be the vector of unknown coefficients at other locations. Then Ac ≤ 0 can be decomposed as

Ac = F 0 G G ′ 0 H • c * c ★ ≤ 0
where Fc * ≤ 0 include the inequalities that involve only variables in c * , Gc * + G ′ c ★ ≤ 0 include the inequalities that involve both c * and c ★ , and Hc ★ ≤ 0 include the inequalities that involve only variables in c ★ . Since Hc ★ ≤ 0 is not related to c * , our first (slight) improvement is to consider only

Bc * ≤ b with B := F G and b := 0 -G ′ • c ★ (2)
in the projection onto c ★ . Furthermore, in the projection cone proj(B):

y * y ★ ≥ 0, y T * , y T ★ • F G = 0, (3)
we only need to consider the generators for (y T * , y T ★ ) T projected onto the dimensions y ★ that correspond to the rows of G since the vector b in (2) has all zero at the dimensions y * . Thus, our second (major) improvement is to project the projection cone proj(B) onto y ★ , and then compute the generators of the projected cone on the dimensions y ★ instead of on the full dimension (y T * , y T ★ ) T .

Correctness and Accuracy

We show that the improvements through our location-by-location idea have the same accuracy as compared with the original work [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], thus generate correct linear invariants. This can be easily observed as follows. Suppose that the CNF formula obtained from Step A4 is Φ = 𝑖 𝑁 𝑖 𝑗=1 𝑐 𝑖,𝑗 and we denote by ⟦Φ⟧ the set of all real vectors (as valuations over all unknown coefficients in the template) that satisfy Φ. The approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] computes a set G of generators by first expanding the CNF Φ into DNF and then collecting the generators of every unsubsumped disjunctive clause of the DNF. Hence, from [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] we have that the convex cone generated by G is equal to that of ⟦Φ⟧. In our improvements, we have that (i) the expansion reordering does not lose accuracy since the reordered CNF Φ ′ is logically equivalent to Φ, and (ii) the projection based generator computation preserves accuracy since always the same generators will be produced by either (i) first computing the generators then performing the projection onto the target location or (ii) first performing the projection then computing the generators. Hence the convex cone generated by the set of generators computed by our improvements is equal to that by ⟦Φ⟧ projected to the target location.

SEGMENTED SUBSUMPTION TESTING

In this part, we propose an independent improvement on the CNF-to-DNF expansion in Step A4 of Section 3.1. Recall that in Section 4.1, we have already proposed an improvement on the CNF-to-DNF by reordering of the expansion order to detect possible subsumptions earlier. The improvement proposed in this part, called segmented subsumption, follows a different idea: we consider to detect local subsumptions in a segment of the whole CNF, so that such local subsumptions can be detected in advance before taking the improvements described in Section 4.

To illustrate our improvement, we depict the general CNF-to-DNF expansion in Figure 2. In the figure, we represent the original CNF formula as 𝑚 𝑖=1 𝐶 𝑖 where we have that each 𝐶 𝑖 = 𝑁 𝑖 𝑗=1 𝑐 𝑖,𝑗 , with each 𝑐 𝑖,𝑗 being a polyhedral cone; we render the conjunctive clauses 𝐶 1 , . . . , 𝐶 𝑚 in the row-wise fashion, for which the vertical order from 𝐶 1 to 𝐶 𝑚 is the order in the CNF-to-DNF expansion and the arrows represent the directions of the expansion in one step. We also define the notion of paths based on Figure 2. A path 𝜋 in Figure 2 is of the form

𝑐 𝑝,𝑘 𝑝 → • • • → 𝑐 𝑖,𝑘 𝑖 → • • • → 𝑐 𝑞,𝑘 𝑞
where 1 ≤ 𝑝 ≤ 𝑞 ≤ 𝑚 and 1 ≤ 𝑘 𝑖 ≤ 𝑁 𝑖 for all 𝑝 ≤ 𝑖 ≤ 𝑞. In other words, a path reflects a partial expansion from 𝐶 𝑝 to 𝐶 𝑞 where the literals chosen in the expansion are 𝑐 𝑝,𝑘 𝑝 , . . . , 𝑐 𝑖,𝑘 𝑖 , . . . , 𝑐 𝑞,𝑘 𝑞 . Note that if 𝑝 = 1 and 𝑞 = 𝑚, then we call the path 𝜋 complete and it corresponds to a disjunctive clause ∧ 𝑞 𝑖=𝑝 𝑐 𝑖,𝑘 𝑖 in the transformed DNF, otherwise it is incomplete. Below we describe the main idea of our segmented subsumption by an illustrative example.

Example 5. Consider a specific scenario of Figure 2 that 𝑚 = 4 and 𝑘 𝑖 = 2 for all 1 ≤ 𝑖 ≤ 4. This scenario is depicted in Figure 3. In this scenario, we consider the situation that the complete path 𝜋 = 𝑐 1,1 → 𝑐 2,1 → 𝑐 3,1 → 𝑐 4,1 subsumes all other complete paths in the expansion, i.e., 𝑖 𝑐 𝑖,𝑘 𝑖 ⊆ 𝑖 𝑐 𝑖,1 for all the choices of 𝑘 𝑖 's. In the original subsumption testing from [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], the other complete paths are expanded and are checked subsumption with the first complete path 𝜋 in the expansion; thus, there could be 15 subsumptions detected by the subsumption testing. Note the reordering in Section 4.1 is irrelevant here as we can assume that the order 𝐶 1 to 𝐶 4 in Figure 3 is already reordered. Although the original approach in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] successfully detected all the subsumptions, it would still perform 15 subsumption testings, which would be time-consuming when there is a large amount of such subsumptions.

In order to further reduce the number of subsumption testing as shown above, we consider to detect local subsumptions in advance. In detail, for the scenario of Figure 3, our local subsumption testing is as follows. We first divide 𝐶 1 , 𝐶 2 , 𝐶 3 , 𝐶 4 into two segmments, 𝐶 1 ∧ 𝐶 2 and 𝐶 3 ∧ 𝐶 4 . Then in each segment, we perform local subsumption testing. For example, for the segment {𝐶 1 , 𝐶 2 }, the local subsumption testing first expand 𝑐 1,1 ∧ 𝑐 2,1 and then checks whether 𝑖=3 𝑐 𝑖,𝑘 𝑖 ⊆ 𝑐 3,1 ∧𝑐 4,1 by 3 subsumption testing, so that only the incomplete path 𝑐 3,1 → 𝑐 4,1 is kept. Thus in this scenario, local subsumtion testing with two segments (𝐶 1 ∧ 𝐶 2 and 𝐶 3 ∧ 𝐶 4 ) may perform only 6 subsumption testing instead of 15 subsumption testing, hence has the potential to reduce the amount of runtime by reducing the number of possible subsumption testing. □ In summary, we divide the whole conjunctive clauses 𝐶 1 , . . . , 𝐶 𝑚 into a number of segments, and perform local subsumption testing to detect subsumption within each segment in advance, so that the overall number of subsumption testing may be improved. Note that different from the original subsumption testing in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that converts a polyhedral cone 𝑃 = 𝑖 𝑐 𝑖,𝑘 𝑖 into the generated invariants 𝐼 by taking the convex closure of 𝑃 and 𝐼 during the CNF-to-DNF expansion, the local subsumption testing in each segment cannot take the covex closure as it may cause inaccurate over-approximation to the original CNF. Thus, to keep the accuracy, in local subsumption testing we only have pairwise subsumption testing between different conjunctions expanded from the local segment (e.g., in the scenario of Figure 3 we have pairwise subsumption testing between 𝑐 1,1 ∧ 𝑐 2,1 , 𝑐 1,1 ∧ 𝑐 2,2 , 𝑐 1,2 ∧ 𝑐 2,1 , 𝑐 1,2 ∧ 𝑐 2,2 within the segment 𝐶 1 ∧ 𝐶 2 ). In practice, we keep the size of each segment small to avoid the combinatorial explosion from the expansion of each segment. clauses could also be segmented with different sizes, but for simplicity we choose a equal size (except for the last segment).

Step B2 In the second step, we perform local subsumption testing at every segment. Given a segment ), then we remove the subsumed one (e.g., 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ) from the DNF. After this step, we get a DNF obtained from ( 4) by removing the subsumed 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 's.

𝐶 𝑝 ∧ • • • ∧ 𝐶 𝑞 with 𝐶 𝑖 = 𝑐 𝑖,1 ∨ • • • ∨ 𝑐 𝑖,
Step B3 In the last step, we group the DNFs obtained from Step B2 for all segments conjunctively together to form a new CNF Φ ′ . (Note that this new CNF Φ ′ is logically equivalent to the original CNF Φ in Step B1.) We then run our invariant-generation procedure from the previous section (Section 4.2) (that incorporates various improvements related to the idea of generating the invariants one program location at a time) over the input Φ ′ to derive the invariants. Since Φ is logically equivalent to Φ ′ , the obtained invariants are guaranteed to be correct.

The integration of the subsumption testing with our location-by-location idea is to first perform the reordered expansion in Section 4.1, then segmented subsumption testing, and finally apply the CNF-to-DNF expansion with projection-based generator computation from Section 4.2 to the CNF formula after segmented subsumption. The pseudo-code for our segmented subsumption testing is given in Algorithm 2 of Appendix A.

EXPERIMENTAL RESULTS

Implementation. We implemented our algorithms as an extension of StInG [START_REF] Sting | StInG: Stanford Invariant Generator[END_REF] in C++, and used PPL 1.2 [START_REF] Bagnara | Possibly Not Closed Convex Polyhedra and the Parma Polyhedra Library[END_REF] for polyhedra manipulation (e.g., projection, generation computation, etc.), and our algorithms work on the original input format of StInG [START_REF] Sting | StInG: Stanford Invariant Generator[END_REF]. To have a prototype implementation that synthesizes linear invariants directly over programs, we wrote most of the benchmark examples in C programming language by having statements that specify the change of values on each transitions, and implemented the transformation from these examples into the format in StInG [START_REF] Sting | StInG: Stanford Invariant Generator[END_REF] via Sparse 0.6.4 [START_REF]Sparse: C language semantic parser[END_REF] (a C language semantic parser). An example input for one of our benchmarks (Scheduler) written in C is given in Appendix C. All experimental results were obtained on an Intel Core i7-7700 (3.6 GHz) machine with 15.5 GiB of memory, running Ubuntu 20.04.2 LTS. Benchmarks. We consider benchmarks from a variety of application domains [START_REF]Arduino: An open-source electronics platform based on easy-to-use hardware and software[END_REF][START_REF] Lamport | A Fast Mutual Exclusion Algorithm[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF][START_REF] Sting | StInG: Stanford Invariant Generator[END_REF][START_REF] Sun | Verifying Stateful Timed CSP Using Implicit Clocks and Zone Abstraction[END_REF]:

Scheduler. The invariant analysis for cooperative multitasking (task scheduling) activated by interrupts and pre-emptive programming (in e.g. Arduino [START_REF]Arduino: An open-source electronics platform based on easy-to-use hardware and software[END_REF]) can be used to ensure the liveness of scheduling. We consider related benchmarks taken from [START_REF] Halbwachs | Verification of Real-Time Systems using Linear Relation Analysis[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. The benchmarks here are divided into two categories, namely "Scheduler'" and "Scheduler". The category "Scheduler'" contains original benchmarks from [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] (with original input format to StInG) that follow a non-standard setting to assign an initial condition to every location in the LinTS in order to model the scenario that the initial location is nondeterministic. The category "Scheduler" considers these benchmarks from StInG under the standard setting of a fixed initial location, corrects several places in the benchmarks that deviate from the original description in [START_REF] Halbwachs | Verification of Real-Time Systems using Linear Relation Analysis[END_REF], and writes them in C.

Fischer. The Fischer mutual exclusion protocol [START_REF] Lamport | A Fast Mutual Exclusion Algorithm[END_REF] is a timing-based algorithm that implements mutual exclusion for a distributed system with skewed clocks, and the invariant analysis can be used to adjust some parameters to ensure critical reachability properties. We take related benchmarks from the PAT toolkit [START_REF] Sun | Verifying Stateful Timed CSP Using Implicit Clocks and Zone Abstraction[END_REF] and rewrite them in C. Cars. A scenario of car systems is illustrated as a dynamic decision problem [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] in which invariant analysis can be used to ensure the safety of autopilot. A car system has 𝑛 cars on a straight road and their acceleration and velocity are determining by their controllers. The lead car starts at an arbitary acceleration with initial velocity 0. The controllers of other cars detect whether the distance between the front car and itself is too close or too far, and adjust their acceleration accordingly. We take the benchmarks from [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], modify them more naturally as LinTS with multiple locations (that distinguish different status on the distance between adjacent cars), and write them in C. We name the category of these benchmarks as "Cars".

For each benchmark above, we consider a variety on the number 𝑟 of processes denoted by "𝑟 -p" (e.g. "2p" stands for two processes) in the benchmark. For all the benchmarks, we compare the running time between our algorithms and the original algorithm in StInG. Experimental Results. Our experimental results are summarized in Table 3 -Table 6. Due to the limit of space, we only present the most representative experimental results and relegate others in Appendix B. Below we first describe the tables, then discuss the experimental results in detail. Table Description. The detailed description of the tables is as follows. In all the tables, for the time consumption and the various size quantities of each benchmark, we have that "Loc" means the number of locations, "Dim" means the number of total dimensions (i.e., the number of all unknown coefficients at all locations), "Line" means the number of lines of the code representation in C, "Time" is the total amount of runtime (by summing up the time consumption at all locations) measured in seconds, "Banged" means the number of successful subsumptions (in the CNF-to-DNF expansion), "Our Approach" means the experimental results by our approach, "StInG" means the experimental results by StInG (which we adapt to PPL 1.2 [START_REF] Bagnara | Possibly Not Closed Convex Polyhedra and the Parma Polyhedra Library[END_REF]), the symbol "-" means not applicable due to either the absence of code representation or time-out or out-of-memory, and "Speedup" in each table shows the ratio of the amount of runtime consumed by StInG against our corresponding approach specified in the table. We set a time bound of 10 hours for time-out in all the tables. Table 3 presents the experimental results on our reordered CNF-to-DNF expansion (Section 4.1). We evaluate the expansion order ( †) that expands (i) intra-transitions from the target location to itself first, then (ii) inter-transitions that involve the target location and one non-target location second, next (iii) the inter-transitions between different non-target locations third, and finally (iv) the intra-transitions over each individual non-target locations last. The expansion order follows the heuristics that the intra-transitions over the target location are the most relevant, and the counterpart over non-target locations are the least relevant. For the projection needed to produce the generators at the target location, we choose the basic method "Gen-Proj" that corresponds to the paragraph Generator Projection in Section 4.2.

Table 4 presents the experimental results on projection-based generator computation (Section 4.2) for which we use the expansion order in Table 3. Furthermore, the tags "Gen-Proj", "FME", "Block", "Block+" correspond to the methods from the paragraphs Generator Projection, Fourier-Motzkin The tag "Gen-Time" records the amount of time (in seconds) consumed on generator computation. Table 5 presents the experimental results on segmented subsumption testing (Section 5) with our reordered expansion in Table 3 and "Block+" in Table 4. In the table, we use the tag "none" to indicate the absence of segmented subsumption, and the tags "Two" (resp. "Three", "Four") to indicate that the size of each segment is two (resp. three, four) except for the last segment, respectively.

Table 6 presents experimental results on parallel computation bestowed from our location-bylocation idea. For the parallel computation, we assume the scenario that every location is assigned a distinct processor to generate the invariants at the location, so that the overall runtime here is the maximal runtime among all locations. In the table, we evaluate our approach with reordered expansion from Table 3 and "Block+" from Table 4, under the setting of both without (indicated by "none") and with the segmented subsumption testing of segment size two (indicated by "Two"); the column "Max" records the maximal runtime among all locations, the "Speedup" column shows the ratio of the amount of runtime by StInG against the "Max" column, and the "Speedup(original)" column records the ratio of the amount of runtime by StInG against the summation of runtime at at all locations (i.e., without parallelism).

Discussion. In all the benchmarks, the generated invariants by our approach (in each table) are exactly the same as from StInG, thus we focus on the comparison in runtime.

Table 3 compares the runtime under our reordered CNF-to-DNF expansion with the most basic projection method of "Gen-Proj" and without segmented subsumption testing. From the table, one can observe that our reordered expansion improves the performance of StInG up to orders of magnitude. (Note that the speedup value with > means that StInG times out on the benchmark). This demonstrates that our reordered expansion can indeed improve the scalability. The improvement can also be observed by the fact that the expansion reordering induces a smaller number of successful subsumptions, which means that these subsumptions are detected earlier by the reordering.

Based on our reordered CNF-to-DNF expansion, Table 4 further compares the runtime from various projection-based generator-computation methods. From the table, one can observe that the simplest way "Gen-Proj" suffices to have comparable speedup against other methods when the dimension is moderate. However, for several high dimensional examples (such as Scheduler-7p and Scheduler-7p'), "Block" and "Block+" are the only methods that help complete the invariant generation in tractable time. Moreover, when comparing "Block" and "Block+", one observes that "Block+" consistently improves "Block" in the runtime for generator computation. This shows that our "Block+" improvement can further leverage the scalability for high dimensional examples in addition to reordered expansion.

Table 5 shows that in addition to our reordered CNF-to-DNF expansion ( †) and "Block+", segmented subsumption testing is effective in high dimensional examples, such as Scheduler-7p, Scheduler'-7p and from Fischer-8p to Fischer-13p, especially for the Fischer benchmarks where a large amount of local subsumption exists, segmented subsumption testing makes it possible to scale the invariant generation to Fischer-12p and Fischer-13p. However, incorporating segmented subsumption is not always better since itself incurs extra runtime at local subsumption testing.

Finally, Table 6 checks the advantage of parallel processing endowed by our location-by-location idea implemented with our reordered CNF-to-DNF expansion ( †) and "Block+". Under the setting of with and without segmented subsumption. From the experimental data, one observes that parallel processing further increases the scalability by a factor between 2 and 3.

Note that our approach cannot scale beyond Cars-4p since we encounter out-of-memory, although the runtime for Cars-2p, Cars-3p and Cars-4p is relatively low. Since our implementation is an extension of StInG, we believe that by improving the memory management in StInG, our approach could handle larger examples from this class of benchmarks.

Remark 2. Here we compare our experiment results with abstract interpretation and recurrence analysis. On one hand, in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], abstract interpretation with standard polyhedral abstract domain [START_REF] Bagnara | Precise Widening Operators for Convex Polyhedra[END_REF][START_REF] Cousot | Automatic Discovery of Linear Restraints Among Variables of a Program[END_REF] has been compared with the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that implements Farkas' Lemma with several reasonale heuristics; the comparison was that the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] (with heuristics that incurs accuracy loss) often generates much more accurate invariants than abstract interpretation and can handle examples where abstract interpretation seemingly diverges. Note that recent advances [START_REF] Singh | Fast polyhedra abstract domain[END_REF] in polyhedral abstract domain explores possible speed up through the separation of independent variables (which is not the case in our benchmarks since all variables in our benchmarks are correlated) , but do not improve the accuracy [START_REF] Singh | Fast polyhedra abstract domain[END_REF]. On the other hand, the state of the art recurrence solver OCRS [START_REF] Ocrs | OCRS: Operational calculus recurrence solver[END_REF] can only generate invariants that relates the value of a program variable and the loop counter, thus cannot be applied to our benchmarks. □ Remark 3. We remark that popular software verification frameworks such as CPAchecker [START_REF] Cpachecker | CPAchecker: The Configurable Software-Verification Platform[END_REF], SeaHorn [START_REF] Seahorn | SeaHorn: A fully automated analysis framework for LLVM-based languages[END_REF] and Ultimate Automizer [START_REF]UltimateAutomizer: A Software Model Checker[END_REF] only consider to check the validity of a given assertion at a program location, thus invariant generation in these frameworks requires an additional assertion as the post-condition and only aims at generating enough invariants to prove/disprove the assertion. In contrast, we follow the setting in [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Cousot | Automatic Discovery of Linear Restraints Among Variables of a Program[END_REF][START_REF] Halbwachs | Verification of Real-Time Systems using Linear Relation Analysis[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that does not require post-condition and aims at generating as much invariants as possible. Thus, these frameworks do not apply to our case. Nevertheless, our improvements can still be applied to the case with assertions, as one can focus on the target program location where the assertion lies, and perform segmented subsumption to simplify the CNF formula. □

RELATED WORKS

Below we compare our approach with most related approaches on numerical invariant generation. We present the comparison classifying the related approaches into different categories of methods. Constraint Solving. Our approach falls in this category. Since we focus on linear invariant generation, our approach is incomparable with those for polynomial invariant generation [START_REF] Adjé | Property-based Polynomial Invariant Generation Using Sums-of-Squares Optimization[END_REF][START_REF] Chatterjee | Polynomial invariant generation for non-deterministic recursive programs[END_REF][START_REF] Chen | Counterexample-Guided Polynomial Loop Invariant Generation by Lagrange Interpolation[END_REF][START_REF] Cousot | Proving Program Invariance and Termination by Parametric Abstraction, Lagrangian Relaxation and Semidefinite Programming[END_REF][START_REF] Steven De Oliveira | Polynomial Invariants by Linear Algebra[END_REF][START_REF] Hrushovski | Polynomial Invariants for Affine Programs[END_REF][START_REF] Humenberger | Automated Generation of Non-Linear Loop Invariants Utilizing Hypergeometric Sequences[END_REF][START_REF] Kapur | Automatically Generating Loop Invariants Using Quantifier Elimination[END_REF][START_REF] Lin | Proving total correctness and generating preconditions for loop programs via symbolic-numeric computation methods[END_REF][START_REF] Rodríguez | Automatic Generation of Polynomial Loop Invariants: Algebraic Foundations[END_REF][START_REF] Sriram Sankaranarayanan | Non-linear loop invariant generation using Gröbner bases[END_REF][START_REF] Yang | Recent advances in program verification through computer algebra[END_REF]. Below we compare our approach with the approaches for linear invariant generation in constraint solving [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF][START_REF] Steven De Oliveira | Synthesizing Invariants by Solving Solvable Loops[END_REF][START_REF] Gupta | InvGen: An Efficient Invariant Generator[END_REF][START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]. The approach [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF] first establishes the framework of linear invariant generation through Farkas' Lemma and solves the invariant through the complete method of quantifier elimination. Quantifier elimination usually has high runtime complexity and is impractical even for programs in moderate size [START_REF] Zhao | A "Hybrid" Approach for Synthesizing Optimal Controllers of Hybrid Systems: A Case Study of the Oil Pump Industrial Example[END_REF]. Thus, our approach has much better runtime performance by developing various techniques to improve scalability.

The approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] considers several heuristics to solve the non-linear constraints from the application of Farkas' Lemma, hence is more scalable than [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF]. A main disadvantage of [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] is that it solves the invariants at all program locations in a single invariant-generation process, thus causing a potentially large amount of combinatorial explosion. Our approach is based on [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] and further incoporates (i) the novel idea of handling the program locations one by one with reinforcing technical improvements and (ii) the segmented subsumption testing to further discover subsumption in the CNF-to-DNF expansion, thus can substantially mitigate the combinatorial explosion from [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] and achieve much better scalability.

The approach [START_REF] Steven De Oliveira | Synthesizing Invariants by Solving Solvable Loops[END_REF] uses eigenvectors to handle several restricted classes of linear invariants. Our approach tackles the general class of affine programs and invariants through completely different techniques, thus focuses on a completely different aspect. The tool INVGEN [START_REF] Gupta | InvGen: An Efficient Invariant Generator[END_REF] focuses on how one integrates abstract interpretation and constraint solving, while our approach aims at improving scalability of constraint solving. Thus the focuses of our approach and INVGEN are orthogonal. Abstract Interpretation. The most classical method to find inductive invariants is abstract interpretation [START_REF] Adjé | Coupling policy iteration with semi-definite relaxation to compute accurate numerical invariants in static analysis[END_REF][START_REF] Bagnara | Generation of Basic Semi-algebraic Invariants Using Convex Polyhedra[END_REF][START_REF] Chakarov | Expectation Invariants for Probabilistic Program Loops as Fixed Points[END_REF][START_REF] Cousot | The ASTREÉ Analyzer[END_REF][START_REF] Müller | Computing polynomial program invariants[END_REF][START_REF] Rodríguez | An Abstract Interpretation Approach for Automatic Generation of Polynomial Invariants[END_REF][START_REF] Rodríguez | Automatic generation of polynomial invariants of bounded degree using abstract interpretation[END_REF]]. An abstraction-interpretation based approach follows the paradigm of first having an abstract domain for the desired invariants, and then perform forward propagation with widening to reach a fixed point. Compared with constraint solving, abstract interpretation does not provide guarantee on the accuracy of the generated invaraints, except for some rare cases [START_REF] Giacobazzi | Completeness in Abstract Interpretation: A Domain Perspective[END_REF]. This point is experimental observed in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], where the method of constraint solving (even with some reasonable heuristics that cause accuracy loss) can produce much tighter invariants than an abstract interpretation through polyhedral abstract domain. We prove that our approach generates the same invariants as [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], hence inherits the advantage of high accuracy from constraint solving. Recurrence Analysis. Recurrence analysis [START_REF] Breck | Templates and recurrences: better together[END_REF][START_REF] Farzan | Compositional Recurrence Analysis[END_REF][START_REF] Humenberger | Algebra-Based Synthesis of Loops and Their Invariants[END_REF][START_REF] Kincaid | Compositional recurrence analysis revisited[END_REF][START_REF] Kincaid | Non-linear reasoning for invariant synthesis[END_REF] first transforms the invariantgeneration problem into a recurrence relation, and then solves the invariants through the analysis of the recurrence relation. Compared with constraint solving, recurrence-based approaches are only applicable to the situations where the underlying recurrence relation has a closed form solution, while constraint solving can be applied when a closed form solution does not exist. Other Methods. Invariant could also be generated through logical inference [START_REF] Dillig | Inductive invariant generation via abductive inference[END_REF][START_REF] Gan | Nonlinear Craig Interpolant Generation[END_REF][START_REF] Garg | ICE: A Robust Framework for Learning Invariants[END_REF][START_REF] Gulwani | Constraint-Based Invariant Inference over Predicate Abstraction[END_REF][START_REF] Kenneth | Quantified Invariant Generation Using an Interpolating Saturation Prover[END_REF][START_REF] Sharma | From invariant checking to invariant inference using randomized search[END_REF][START_REF] Xu | Interval counterexamples for loop invariant learning[END_REF], machine learning [START_REF] Garg | Learning invariants using decision trees and implication counterexamples[END_REF][START_REF] He | Learning fast and precise numerical analysis[END_REF][START_REF] Yao | Learning nonlinear loop invariants with gated continuous logic networks[END_REF] and dynamic analysis [START_REF] Csallner | DySy: dynamic symbolic execution for invariant inference[END_REF][START_REF] Nguyen | Using dynamic analysis to discover polynomial and array invariants[END_REF][START_REF] Sharma | A Data Driven Approach for Algebraic Loop Invariants[END_REF]. Since these methods either utilize specilized inference rules, or depend on a potentially big dataset, or require a possible large amount of program executions, they could not ensure any accuracy guarantee on the generated invariants.

CONCLUSION AND FUTURE WORK

We proposed two improvements to the previous approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that solves the invariant generation problem with Farkas' Lemma and several reasonable heuristics. The first is a novel idea that generates the linear invariants for each program location separately. The idea enables one to speed up the invariant-generation process by parallel processing, and we reinforce the idea with two technical improvments that further reduce the combinatorial explosion from [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF]: (i) in the CNF-to-DNF expansion, we reordered the expansion so that the target location (over which the invariants are to be generated) comes first to detect subsumptions earlier; (ii) in the generator computation, we proposed various methods that calculate the generators projected onto the target location, including a novel method that improves Block Elimination. The second is a segmented subsumption testing in addition to the subsumption testing in [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF] that could detect subsumptions at local conjunctive clauses. Experimental results show that our improvements can improve the scalability of the approach [START_REF] Sankaranarayanan | Constraint-Based Linear-Relations Analysis[END_REF], the only known practical algorithm that is based on Farkas' Lemma, by orders of magnitude. A future direction would be to investigate new approaches to further improve scalability. Another future direction would be to extend our improvements to polynomial invariant generation [START_REF] Chatterjee | Polynomial invariant generation for non-deterministic recursive programs[END_REF]. A further direction is to integrate the approach in [START_REF] Sharma | Simplifying Loop Invariant Generation Using Splitter Predicates[END_REF] to handle disjunctive invariants.

A PSEUDO-CODE FOR SEGMENTED SUBSUMPTION TESTING Algorithm 2 Segmented Subsumption Testing Input: 𝑚 𝑖=1 𝐶 𝑖 : the CNF formula obtained from Step A4 (Section 3.1) in the approach [START_REF] Michael Colón | Linear Invariant Generation Using Non-linear Constraint Solving[END_REF] where each 𝐶 𝑖 = 𝑁 𝑖 𝑗=1 𝑐 𝑖,𝑗 is the 𝑖-th conjunctive clause and is a disjunction of polyhedra over the unknown coefficients in the template; 𝑘 : the size of segments (the size could be different but for simplicity we choose a equal size except for the last segment) Output: Φ ′ : the segmented CNF formula 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ← expand 𝜙 𝑖 fully into its equivalent DNF; //𝑄 𝑘 𝑝 ,.. Φ 𝑖 ← 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ;

Fig. 1 .

 1 Fig. 1. The LinTS for a Vagrant Robot

, Vol. 1 ,

 1 No. 1, Article . Publication date: March 2022.

Fig. 2 .

 2 Fig. 2. CNF-to-DNF Expansion Fig. 3. Example of Expansion

, Vol. 1 ,

 1 No. 1, Article . Publication date: March 2022.

1 :

 1 {𝜙 0 , • • • , 𝜙 ( ⌈𝑚/𝑘 ⌉-1) } ← divide 𝑚 𝑖=1 𝐶 𝑖 into ⌈𝑚/𝑘⌉ segments; //𝜙 𝑖 = 𝐶 𝑝 ∧ • • • ∧ 𝐶 𝑞 where 𝑝 = 𝑖 • 𝑘 + 1 and 𝑞 = 𝑖 • 𝑘 + 𝑘 ( except 𝑞 = 𝑚 if 𝑖 = ⌈𝑚/𝑘⌉ -1 ) 2:Φ 𝑖 ← ∅ ( for 0 ≤ 𝑖 ≤ ⌈𝑚/𝑘⌉ -1 ); 3: for 𝑖 = 0 to ⌈𝑚/𝑘⌉ -1 do 4:

11: end for 12 :

 12 Φ ′ ← Φ 0 • • • Φ ( ⌈𝑚/𝑘 ⌉-1) ; 13: return Φ ′ ; C AN EXAMPLE OF BENCHMARKS # define time_cnt __attribute__ (( time_cnt )) # define p1_p1 __attribute__ (( p1_p1 )) # define p2_p2 __attribute__ (( p2_p2 )) # define p3_p3 __attribute__ (( p3_p3 )) # define p4_p4 __attribute__ (( p4_p4 )) # define p5_p5 __attribute__ (( p5_p5 )) # define p6_p6 __attribute__ (( p6_p6 )) # define p7_p7 __attribute__ (( p7_p7 )) # define p0_p2 __attribute__ (( p0_p2 )) # define p0_p1 __attribute__ (( p0_p1 )) # define p1_p2 __attribute__ (( p1_p2 )) # define p2_p1 __attribute__ (( p2_p1 )) # define p2_p0 __attribute__ (( p2_p0 )) # define p1_p0 __attribute__ (( p1_p0 )) # define p1_p3 __attribute__ (( p1_p3 )) # define p3_p1 __attribute__ (( p3_p1 )) # define p1_p4 __attribute__ (( p1_p4 )) # define p4_p1 __attribute__ (( p4_p1 )) # define p1_p5 __attribute__ (( p1_p5 )) # define p5_p1 __attribute__ (( p5_p1 )) # define p1_p6 __attribute__ (( p1_p6 )) # define p6_p1 __attribute__ (( p6_p1 )) # define p1_p7 __attribute__ (( p1_p7 )) # define p7_p1 __attribute__ (( p7_p1 )) # define p2_p3 __attribute__ (( p2_p3 )) # define p3_p2 __attribute__ (( p3_p2 )) # define p2_p4 __attribute__ (( p2_p4 )) # define p4_p2 __attribute__ (( p4_p2 )) # define p2_p5 __attribute__ (( p2_p5 )) # define p5_p2 __attribute__ (( p5_p2 )) # define p2_p6 __attribute__ (( p2_p6 )) # define p6_p2 __attribute__ (( p6_p2 )) # define p2_p7 __attribute__ (( p2_p7 )) # define p7_p2 __attribute__ (( p7_p2 )) # define p3_p4 __attribute__ (( p3_p4 )) # define p4_p3 __attribute__ (( p4_p3 )) # define p3_p5 __attribute__ (( p3_p5 )) # define p5_p3 __attribute__ (( p5_p3 )) # define p3_p6 __attribute__ (( p3_p6 )) # define p6_p3 __attribute__ (( p6_p3 )) # define p3_p7 __attribute__ (( p3_p7 )) # define p7_p3 __attribute__ (( p7_p3 ))

  e., 𝜑 |= 𝜓 ) if and only if there exist non-negative real numbers 𝜆 0 , 𝜆 1 , . . . , 𝜆 𝑚 such that (i) 𝑐 𝑗 = 𝑚 𝑖=1 𝜆 𝑖 • 𝑎 𝑖 𝑗 for all 1 ≤ 𝑗 ≤ 𝑛, and (ii) 𝑑 = 𝜆 0 + 𝑚 𝑖=1 𝜆 𝑖 • 𝑏 𝑖 . Moreover, 𝜑 is unsatisfiable if and only if the inequality -1 ≥ 0 (as 𝜓 ) can be derived from above.

, Vol. 1, No. 1, Article . Publication date: March 2022.

Table 1 .

 1 The Tabular Form for Farkas' Lemma

  2 𝑖=1 𝑐 𝑖,𝑘 𝑖 ⊆ 𝑐 1,1 ∧ 𝑐 2,1 for other choices of 𝑘 𝑖 's; if indeed 2 𝑖=1 𝑐 𝑖,𝑘 𝑖 ⊆ 𝑐 1,1 ∧ 𝑐 2,1for other choices of 𝑘 𝑖 's, then only 3 subsumption testing is performed and we only keep the incomplete path 𝑐 1,1 → 𝑐 2,1 . Likewise, for the segment {𝐶 3 , 𝐶 4 }, the local subsumption testing may also detect that4

  𝑁 𝑖 for each 𝑝 ≤ 𝑖 ≤ 𝑞, the detailed local subsumption testing on this segment is as follows. We first expand the segment fully into its equivalent DNF1≤𝑘 𝑝 ≤𝑁 𝑝 ,••• ,1≤𝑘 𝑞 ≤𝑁 𝑞 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 where 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 := Then for each pair of 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 and 𝑄 𝑘 ′ 𝑝 ,...,𝑘 ′ 𝑞 , we check whether one of them (e.g. 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ) is subsumed by the other (e.g. 𝑄 𝑘 ′ 𝑝 ,...,𝑘 ′ 𝑞 ); if the subsumption test is successful (e.g., 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ⊆ 𝑄 𝑘 ′ 𝑝 ,...,𝑘 ′

	𝑞
	𝑐 𝑖,𝑘 𝑖 .	(4)
	𝑖=𝑝
	𝑞

Table 3 .

 3 Experimental Results for the Reordered CNF-to-DNF Expansion ( †)

		Benchmarks				Our Approach
						StInG
							Reordered-Expansion †: Gen-Proj
	Name		Loc Dim Line
						Time ( Banged )	Time ( Banged )	Speedup
		2p	2	16	-	0.01 ( 26 )	0.01 ( 34 )	1.00X
		3p	3	33	-	0.17 ( 380 )	0.10 ( 245 )	1.70X
		4p	4	56	-	60.81 ( 26,629 )	1.38 ( 1,508 )	44.06X
	Scheduler'
		5p	5	85	-	7,436.34 ( 2,548,704 )	26.52 ( 25,996 )	280.40X
		6p	6	120	-	time out	892.75 ( 66,905 )	>40.32X
		7p	7	161	-	time out	time out	-
		3p	3	33	336	0.17 ( 279 )	0.10 ( 261 )	1.70X
		4p	4	56	609	4.16 ( 2,895 )	0.98 ( 1,474 )	4.24X
	Scheduler	5p	5	85	1017	135.80 ( 39,150 )	13.97 ( 8,570 )	9.72X
		6p	6	120 1587	7,541.53 ( 906,454 )	277.96 ( 63,070 )	27.13X
		7p	7	161 2352	time out	time out	-
		6p	7	63	710	9.18 ( 8,423 )	3.07 ( 4,217 )	2.99X
		7p	8	80	987	59.16 ( 32,668 )	13.60 ( 12,566 )	4.35X
		8p	9	99	1327	373.62 ( 127,918 )	70.42 ( 43,191 )	5.30X
	Fischer	9p	10	120 1736	2,345.96 ( 503,369 )	398.26 ( 163,623 )	5.89X
		10p 11	143 2218 14,664.68 ( 1,985,857 )	2,361.87 ( 649,409 )	6.20X
		11p 12	168 2780	time out	14,307.40 ( 2,620,864 )	>2.51X
		12p 13	195 3453	time out	time out	-
		2p	3	27	216	0.01 ( 28 )	0.01 ( 63 )	1.00X
		3p	5	60	616	560.70 ( 788,508 )	1.21 ( 2,299 )	463.38X
	Cars
		4p	7	105 1283	time out	86.83 ( 37,567 )	>414.60X
		5p	9	162 2339	out of memory	out of memory	-

Table 4 .

 4 Experimental Results for Projection-Based Generator Computation

		Benchmarks						Our Approach
						Reordered-Expansion †: Gen-Proj Reordered-Expansion †: FME Reordered-Expansion †: Block Reordered-Expansion †: Block+
	Name		Loc Dim Line
						Time ( Gen-Time )	Speedup	Time ( Gen-Time ) Speedup Time ( Gen-Time ) Speedup Time ( Gen-Time ) Speedup
		2p	2	16	-	0.01 ( 0.01 )	1.00X	0.01 ( 0.01 )	1.00X	0.03 ( 0.02 )	0.33X	0.02 ( 0.02 )	0.50X
		3p	3	33	-	0.10 ( 0.06 )	1.70X	0.14 ( 0.08 )	1.21X	0.31 ( 0.27 )	0.55X	0.14 ( 0.08 )	1.21X
		4p	4	56	-	1.38 ( 0.24 )	44.07X	1.60 ( 0.52 )	38.01X	2.41 ( 1.42 )	25.23X	1.50 ( 0.46 )	40.54X
	Scheduler'	5p	5	85	-	26.52 ( 0.38 )	280.40X	27.11 ( 1.32 )	274.30X	30.49 ( 4.03 )	243.89X	27.27 ( 1.07 )	272.69X
		6p	6	120	-	892.75 ( 38.59 )	>40.32X	906.27 ( 48.88 )	>39.72X	887.27 ( 28.01 )	>40.57X	868.52 ( 11.20 )	>41.45X
		7p	7	161	-	time out ( 61,366.73 )	-	time out ( 66,745.07 )	-	35,017.19 ( 209.25 )	>1.03X	34,960.08 ( 140.49 )	>1.03X
		8p	8	208	-	time out	-	time out	-	time out	-	time out	-
		3p	3	33	336	0.10 ( 0.01 )	1.70X	0.13 ( 0.05 )	1.31X	0.26 ( 0.19 )	0.65X	0.14 ( 0.07 )	1.21X
		4p	4	56	609	0.98 ( 0.17 )	4.24X	1.18 ( 0.39 )	3.53X	1.96 ( 1.14 )	2.12X	1.30 ( 0.37 )	3.20X
		5p	5	85	1017	13.97 ( 1.22 )	9.72X	14.58 ( 1.94 )	9.31X	17.55 ( 4.91 )	7.74X	18.68 ( 1.71 )	7.27X
	Scheduler
		6p	6	120 1587	277.96 ( 35.26 )	27.13X	285.30 ( 45.20 )	26.43X	262.31 ( 21.79 )	28.75X	240.50 ( 6.60 )	31.36X
		7p	7	161 2352 time out ( 61,660.65 )	-	time out ( 66,808.22 )	-	4,759.24 ( 142.23 )	>7.56X	4,698.29 ( 77.68 )	>7.66X
		8p	8	208 3351	time out	-	time out	-	time out	-	time out	-
		6p	7	63	710	3.07 ( 0.09 )	2.99X	3.21 ( 0.40 )	2.86X	3.44 ( 0.59 )	2.67X	3.15 ( 0.21 )	2.91X
		7p	8	80	987	13.60 ( 0.20 )	4.35X	14.18 ( 0.69 )	4.17X	14.73 ( 1.25 )	4.02X	14.91 ( 0.74 )	3.97X
		8p	9	99	1327	70.42 ( 0.35 )	5.31X	71.97 ( 1.47 )	5.19X	73.52 ( 2.66 )	5.08X	71.15 ( 1.26 )	5.25X
	Fischer	9p	10	120 1736	398.26 ( 0.70 )	5.89X	408.79 ( 2.69 )	5.74X	413.76 ( 4.92 )	5.67X	414.08 ( 2.59 )	5.67X
		10p 11	143 2218	2,361.87 ( 1.34 )	6.21X	2,429.26 ( 5.03 )	6.04X	2,401.19 ( 8.75 )	6.11X	2,398.38 ( 4.54 )	6.11X
		11p 12	168 2780	14,307.40 ( 2.25 )	>2.52X	14,469.76 ( 8.58 )	>2.49X	14,483.90 ( 15.02 )	>2.49X	14,318.84 ( 7.96 )	>2.51X
		12p 13	195 3453	time out	-	time out	-	time out	-	time out	-
		2p	3	27	216	0.01 ( 0 )	1.00X	0.02 ( 0 )	0.50X	0.02 ( 0 )	0.50X	0.01 ( 0 )	1.00X
		3p	5	60	616	1.21 ( 0.02 )	463.39X	1.25 ( 0.01 )	448.56X	1.55 ( 0.09 )	361.74X	1.21 ( 0.04 )	463.39X
	Cars
		4p	7	105 1283	86.83 ( 0.11 )	>414.60X	86.43 ( 0.58 )	>416.52X	87.20 ( 1.29 )	>412.84X	84.49 ( 0.56 )	>426.09X
		5p	9	162 2339	out of memory	-	out of memory	-	out of memory	-	out of memory	-

Table 5 .

 5 Experimental Results for Segmented Subsumption

		Benchmarks					Our Approach( Reordered-Expansion †: Block+ )
						none		Two		Three		Four
	Name		Loc Dim Line
						Time ( Banged )	Speedup	Time ( Banged )	Speedup	Time ( Banged )	Speedup Time ( Banged ) Speedup
		2p	2	16	-	0.02 ( 34 )	0.50X	0.03 ( 26 )	0.33X	0.02 ( 38 )	0.50X	0.02 ( 73 )	0.50X
		3p	3	33	-	0.14 ( 245 )	1.21X	0.20 ( 207 )	0.85X	0.24 ( 195 )	0.71X	0.23 ( 256 )	0.74X
		4p	4	56	-	1.50 ( 1,508 )	40.54X	1.84 ( 1,348 )	33.05X	2.02 ( 1,301 )	30.10X	2.69 ( 1,310 )	22.61X
	Scheduler'	5p	5	85	-	27.27 ( 25,996 )	272.69X	20.86 ( 15,172 )	356.49X	18.69 ( 9,991 )	397.88X	24.35 ( 9,488 )	305.39X
		6p	6	120	-	868.52 ( 66,905 )	>41.45X	720.78 ( 54,862 )	>49.95X	648.79 ( 46,226 )	>55.49X	2,470.41 ( 39,207 )	>14.57X
		7p	7	161	-	34,960.08 ( 629,977 )	>1.03X	19,721.28 ( 437,681 )	>1.83X	26,420.75 ( 449,445 )	>1.36X	time out	-
		8p	8	208	-	time out	-	time out	-	time out	-	time out	-
		3p	3	33	336	0.14 ( 261 )	1.21X	0.19 ( 234 )	0.89X	0.20 ( 212 )	0.85X	0.21 ( 261 )	0.81X
		4p	4	56	609	1.30 ( 1,474 )	3.20X	1.78 ( 1,318 )	2.34X	2.04 ( 1,270 )	2.04X	3.59 ( 1,315 )	1.16X
		5p	5	85	1017	18.68 ( 8,570 )	7.27X	15.21 ( 6,728 )	8.93X	37.73 ( 6,307 )	3.60X	72.60 ( 6,065 )	1.87X
	Scheduler
		6p	6	120 1587	240.50 ( 63,070 )	31.36X	201.53 ( 51,342 )	37.42X	209.50 ( 42,837 )	36.00X	596.82 ( 36,737 )	12.64X
		7p	7	161 2352	4,698.29 ( 492,751 )	>7.66X	2,752.66 ( 299,471 )	>13.08X	3,760.55 ( 314,398 )	>9.57X	8,986.91 ( 263,933 )	>4.01X
		8p	8	208 3351	time out	-	time out	-	time out	-	time out	-
		6p	7	63	710	3.15 ( 4,217 )	2.91X	3.19 ( 2,808 )	2.88X	3.40 ( 2,336 )	2.70X	4.40 ( 2,280 )	2.09X
		7p	8	80	987	14.91 ( 12,566 )	3.97X	9.45 ( 6,188 )	6.26X	8.78 ( 4,407 )	6.74X	10.87 ( 4,270 )	5.44X
		8p	9	99	1327	71.15 ( 43,191 )	5.25X	28.29 ( 14,231 )	13.21X	25.39 ( 10,383 )	14.72X	28.19 ( 8,690 )	13.25X
		9p	10	120 1736	414.08 ( 163,623 )	5.67X	91.96 ( 34,033 )	25.51X	72.05 ( 22,273 )	32.56X	75.95 ( 18,629 )	30.89X
	Fischer	10p 11	143 2218	2,398.38 ( 649,409 )	6.11X	335.47 ( 90,362 )	43.71X	185.65 ( 42,349 )	78.99X	190.35 ( 36,250 )	77.04X
		11p 12	168 2780 14,318.84 ( 2,620,864 )	>2.51X	1,280.06 ( 235,812 )	>28.12X	765.02 ( 129,144 )	>47.06X	583.50 ( 81,593 )	>61.70X
		12p 13	195 3453	time out	-	5,521.77 ( 681,865 )	>6.52X	2,689.64 ( 296,297 )	>13.38X 1,976.83 ( 190,150 ) >18.21X
		13p 14	224 4397	time out	-	23,812.01 ( 1,835,173 )	>1.51X	8,430.01 ( 588,289 )	>4.27X	8,239.66 ( 501,438 )	>4.37X
		14p 15	255 5039	time out	-	time out	-	time out	-	time out	-
		2p	3	27	216	0.01 ( 63 )	1.00X	0.04 ( 61 )	0.25X	0.04 ( 62 )	0.25X	0.05 ( 65 )	0.20X
		3p	5	60	616	1.21 ( 2,299 )	463.39X	2.01 ( 2,279 )	278.96X	1.63 ( 1,592 )	343.99X	2.04 ( 1,612 )	274.85X
	Cars
		4p	7	105 1283	84.49 ( 37,567 )	>426.09X	100.35 ( 37,525 )	>358.74X	89.78 ( 25,423 )	>400.98X	103.59 ( 25,771 )	>347.52X
		5p	9	162 2339	out of memory	-	out of memory	-	out of memory	-	out of memory	-

Table 6 .

 6 Experimental Results for Parallel Computation

		Benchmarks				Our Approach ( Reordered-Expansion †: Block+ )
							none			Two
	Name		Loc Dim Line
						Max	Speedup(original) Speedup	Max	Speedup(original) Speedup
		2p	2	16	-	0.01	0.50X	1.00X	0.02	0.33X	0.50X
		3p	3	33	-	0.05	1.21X	3.40X	0.08	0.85X	2.13X
		4p	4	56	-	0.47	40.54X	129.38X	0.56	33.05X	108.59X
	Scheduler'
		5p	5	85	-	7.52	272.69X	988.87X	7.33	356.49X	1,014.51X
		6p	6	120	-	308.48	>41.45X	>116.70X	302.22	>49.95X	>119.12X
		7p	7	161	-	11,134.30	>1.03X	>3.23X	5,944.38	>1.83X	>6.06X
		3p	3	33	336	0.04	1.21X	4.25X	0.07	0.89X	2.43X
		4p	4	56	609	0.39	3.20X	10.66X	0.49	2.34X	8.49X
	Scheduler	5p	5	85	1017	4.90	7.27X	27.71X	4.42	8.93X	30.72X
		6p	6	120 1587	72.07	31.36X	104.64X	70.04	37.42X	107.67X
		7p	7	161 2352 1,336.13	>7.66X	>26.94X	750.41	>13.08X	>47.97X
		6p	7	63	710	1.22	2.91X	7.52X	0.79	2.88X	11.62X
		7p	8	80	987	7.75	3.97X	7.63X	2.85	6.26X	20.76X
		8p	9	99	1327	45.41	5.25X	8.22X	11.14	13.21X	33.54X
		9p	10	120 1736	283.95	5.67X	8.26X	43.96	25.51X	53.37X
	Fischer
		10p 11	143 2218 1,700.51	6.11X	8.62X	187.85	43.71X	78.07X
		11p 12	168 2780 10,230.51	>2.51X	>3.51X	768.56	>28.12X	>46.84X
		12p 13	195 3453 time out	-	-	3,437.67	>6.52X	>10.47X
		13p 14	224 4397 time out	-	-	14,885.30	>1.51X	>2.42X
			3	27	216	0.01	1.00X	1.00X	0.02	0.25X	0.50X
	Cars	3p	5	60	616	0.83	463.39X	675.54X	1.04	278.96X	539.13X
		4p	7	105 1283	75.83	>426.09X	>474.74X	84.46	>358.74X	>426.24X

Elimination, Block Elimination and Our improvment on Block Elimination in Section 4.2, respectively.

  .,𝑘 𝑞 = 𝑞 𝑖=𝑝 𝑐 𝑖,𝑘 𝑖 where 1 ≤ 𝑘 𝑝 ≤ 𝑛 𝑝 , • • • , 1 ≤ 𝑘 𝑞 ≤ 𝑛 𝑞 5: for each 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 and 𝑄 𝑘 ′ 𝑝 ,...,𝑘 ′ 𝑞 such that 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ≠ 𝑄 𝑘 ′ 𝑝 ,...,𝑘 ′ 𝑞 do //1 ≤ 𝑘 𝑖 , 𝑘 ′ 𝑖 ≤ 𝑁 𝑖 where 𝑝 ≤ 𝑖 ≤ 𝑞 6: if 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ⊆ 𝑄 𝑘 ′ 𝑝 ,...,𝑘 ′ 𝑞 then 7:remove subsumed 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 from 𝑄 𝑘 𝑝 ,...,𝑘 𝑞 ;

	8:	end if
	9:	end for
	10:

, Vol. 1, No. 1, Article . Publication date: March 2022.

B DETAILED EXPERIMENTAL RESULTS