extfieldsofstudy:
Computer Science
Mathematics

s2fieldsofstudy:
Computer Science
Mathematics

text:
A Strongly Exponential Separation of DNNFs from CNF Formulas

Decomposable Negation Normal Forms (DNNFs) are Boolean circuits in negation normal form where the subcircuits leading into each AND gate are defined on disjoint sets of variables. We prove a strongly exponential lower bound on the size of DNNFs for a class of CNF formulas built from expander graphs. As a corollary, we obtain a strongly exponential separation between DNNFs and CNF formulas in prime implicates form. This settles an open problem in the area of knowledge compilation (Darwiche and Marquis, 2002).

Introduction
The aim of knowledge compilation is to succinctly represent propositional knowledge bases in a format that allows for answering a number of queries in polynomial time [DM02]. Choosing a representation language generally involves a trade-off between succinctness and the range of queries that can be efficiently answered. For instance, CNF formulas are more succinct than prime implicate formulas (PIs), but the latter representation enjoys clause entailment checks in polynomial time whereas CNF formulas in general do not, unless P = NP [GKPS95,CM78]. The need to balance the competing requirements of succinctness and tractability has led to the introduction of a large variety of representation languages that strike this balance in different ways.
Decomposable Negation Normal Forms (DNNFs) are a restricted form of Boolean circuits in negation normal form (NNF) such that the subcircuits leading into an AND gate are defined on disjoint sets of variables. DNNFs are among the most succinct representation languages considered in knowledge compilation-for instance, they include variants of binary decision diagrams such as ordered binary decision diagrams (OBDDs) and even free binary decision diagrams (FBDDs, also known as read-once branching programs) as special cases. In spite of their expressive power, they still support clause entailment queries and satisfiability testing in linear time, both of which are of critical importance in practice.
In this paper, we study the relative succinctness of DNNFs and CNF formulas. On the one hand, DNNFs are known to be exponentially more succinct than CNF formulas in some cases [GKPS95]; on the other hand, it was known that CNFs do not have polynomial DNNF representations unless the polynomial hierarchy collapses [SK96,CDLS02]. The main contribution of this paper is in improving this conditional, superpolynomial lower bound to an unconditional, exponential lower bound (Theorem 5): There is a class C of CNFs such that the DNNF size of each F ∈ C is 2 Ω(size(F )) .
We prove this result by combining the combinatorial properties of expander graphs with recently developed techniques for proving lower bounds on the size of non-deterministic read-once branching programs for CNF formulas [Raz14b]. In addition to resolving an open problem by Darwiche and Marquis [DM02], our lower bound has several other implications (cf. Appendix 7): • The formulas used to prove Theorem 5 are PIs, so it immediately follows that PIs can be exponentially more succinct than DNNFs (Corollary 3). Again, this answers an open question by Darwiche and Marquis [DM02].
Related Work. To the best of our knowledge, ours is the first unconditional lower bound for general DNNFs. However, there is a rich literature on lower bounds for more restricted representation languages, such as OBDDs or FBDDs [Weg00,Juk12]. Moreover, certain subclasses of DNNFs, so-called decision-DNNFs, have been recently considered in database theory in the context of probabilistic databases. In this setting, lower bounds are obtained by a quasipolynomial simulation of decision-DNNFs by FBDDs in combination with known exponential lower bounds for FBDDs [BLRS13,BLRS14]. Pipatsrisawat and Darwiche have proposed a framework for showing lower bounds on structured DNNFs [PD10], a subclass of DNNFs in which the variables respect a common tree-ordering.
The formulas used to prove Theorem 5 are based on expander graphs and were originally introduced to establish an exponential lower bound for the OBDD size of CNFs [BS14]. The present paper leverages a recent result by Razgon [Raz14b,Theorem 4] to lift this lower bound to DNNFs; indeed, we even generalize such result, at the same time providing a significantly shorter proof based on a new combinatorial result (Theorem 3) which may be of independent interest.

Preliminaries
For every integer n ≥ 1, we let [n] denote the set {1, . . . , n}. Let X be a countable set of elements, called variables, and call the elements of {0, 1} constants. A literal is a variable in X (a positive literal) or the negation, ¬, of a variable in X (a negative literal). For a literal l ∈ {x, ¬x} for some x ∈ X, we let var(l) = x. An atom is a literal or a constant. An assignment is a function from a subset of X to the constants. We equivalently view an assignment f as {¬x : f (x) = 0} ∪ {x : f (x) = 1}; so, as a set of atoms, an assignment does not contain the constant 0 nor a pair of complementary literals.
Negation Normal Forms (NNFs). A Negation Normal Form (NNF) C on a finite set of variables vars(C) ⊆ X (also known as a De Morgan circuit ) is a node labelled directed acyclic graph (DAG), whose labelled nodes and arcs are respectively called the gates and wires of C. A NNF C has a designated sink (outdegree 0) node, referred to as the output gate of the circuit, and denoted by output(C). The source nodes of C (indegree 0), denoted by inputs(C), are referred to as the input gates of C and are labelled by a constant or a literal in {x, ¬x : x ∈ vars(C)}. The non source nodes of D, referred as internal gates, are labelled by ∧ or ∨. The label of a node v is denoted by label(v), and a node v is also called a label(v)-gate. We let atoms(C) = {label(v) : v ∈ inputs(C)}. In this paper, the size of C, in symbols size(C), is the number of wires in C.
Let G be a DAG and let v be a node in G. The subgraph of G sinked at v is the DAG whose node set is and whose arcs are exactly the arcs of G among the nodes in V ′ . Let C be an NNF, and let v be a node in the DAG G underlying C. We let sub(C, v) denote the subcircuit of C sinked at v, that is, the NNF on variables vars(C) whose underlying DAG is the subgraph of C sinked at v, with the same labels and designated sink v.
Let C be an NNF, v be a gate in C, and f be an assignment of vars(C). The value of sub(C, v) under f , in symbols sub(C, v)(f ), is defined inductively on the structure of sub(C, v) as follows: • if v is a 0-gate or an ∨-gate having fanin 0, then sub(C, v)(f ) = 0; • if v is a 1-gate or a ∧-gate having fanin 0, then sub(C, v)(f ) = 1; The value of C under f , in symbols C(f ), is equal to sub(C, output(C))(f ). We let sat(C) denote the set of all satisfying assignments (or models) of C, that is, assignments f : vars(C) → {0, 1} such that C(f ) = 1. Let C be an NNF on variables vars(C). We say that C computes the (Boolean) function F : ; moreover, C is said trivial if it computes the identically 0 (or 1) function. Two NNFs C and C ′ are equivalent if they compute the same function; in particular, vars(C) = vars(C ′ ).
Let C be an NNF. A certificate of C (sometimes also called proof tree) is an NNF T on variables vars(C) whose gates and wires are subsets of the gates and wires of C satisfying the following: output(C) = output(T ); if a gate v is in C ∩ T , and v is a ∧-gate in C with input wires from gates v 1 , . . . , v i , then the gates v 1 , . . . , v i and the wires (v 1 , v), . . . , and v is an ∨-gate in C with input wires from gates v 1 , . . . , v i , then exactly one gate w ∈ {v 1 , . . . , v i } is in T and the wire (w, v) is in T v (i ≥ 0). We let cert(C) denote the set of certificates of C.
The models and the certificates of an NNF are nicely related as follows.
Proposition 1. Let C be an NNF and let f be an assignment. The following are equivalent.
• There exists T ∈ cert(C) such that f ∈ sat(T ).
Proof. For the forward direction, let f ∈ sat(C). Call a wire (u, v) in the DAG underlying C activated by f if f satisfies the subcircuit of C sinked at u, in symbols sub(C, u)(f ) = 1. It is readily verified that there exists a certificate T for C containing only wires activated by f . Moreover, f ∈ sat(T ) because by construction sub(T, v)(f ) = 1 for all input gates v of T , and therefore sub(T, output(T ))(f ) = 1. For the backward direction, let T be a certificate of C such that f ∈ sat(T ). By induction on the structure of C, we prove that for all t ∈ T it holds that sub(T, t)(f ) = sub(C, t)(f ) = 1. Since output(C) = output(T ), we conclude that C(f ) = sub(C, output(C))(f ) = sub(T, output(T ))(f ) = 1.
If t is an input gate of T , then sub(T, t)(f ) = 1 because otherwise f does not satisfy T . We also have that t is an input gate of C, hence sub(C, t)(f ) = 1. Let t be a ∨-gate in T , with input wires from gates t 1 , . . . , t i in C; say without loss of generality that t 1 is chosen in T . By the induction hypothesis, sub(T, t 1 )(f ) = sub(C, t 1 )(f ) = 1; hence, sub(T, t)(f ) = sub(C, t)(f ) = 1. The case where t is a ∧-gate is similar.
Decomposable Negation Normal Forms (DNNFs). An NNF D is decomposable (in short, a DNNF ) if for all ∧-gates v with input wires from gates v 1 , . . . , v i and all j, j ′ ∈ [i], j = j ′ , the variable sets of the subcircuits of D sinked at v j and v j ′ are disjoint, in symbols, Conjunctive Normal Forms (CNFs). A conjunctive normal form (in short, CNF ) is a finite set of clauses, where a clause is a finite set of literals. Equivalently, we regard a CNF as an NNF where the maximum number of wires on a path from a input gate to the output gate is 2, and each ∨-gate has input wires only from input gates. A CNF F is monotone if its labels do not contain negative literals, a k-CNF if the fanin of ∨-gates is at most k, and a read k times CNF if, for every x ∈ vars(F ), the number of wires leaving nodes whose label contain the variable x is at most k.
Graphs. We refer to a standard reference for basic notions and facts in graph theory [Die05]. Let G = (V, E) be a graph.
For every S ⊆ V , we let N S,G denote the open neighbourhood of S in G, in symbols, (1) Note that a (c, d)-expander is connected; otherwise, if S is the vertex set of a (nonempty) connected component of G of minimum size, then 0 < |S| ≤ |V |/2 and 0 = |N S,G | ≥ c|S| > 0.
Theorem 1 (Section 9.2 in [AS00]). For all d ≥ 3, there exists c > 0 and a sequence of graphs If G has no isolated vertices, we view E as a CNF on the variables vars(E) = V (a trivial CNF, if E = ∅), which we call a graph CNF. Let E be a graph CNF with underlying graph G = (V, E). Note that the satisfying assignments of E correspond to vertex covers of G: If f is a satisfying assignment of E, then {x : f (x) = 1} is a vertex cover of G, and if C is a vertex cover of G, then any assignment f of vars(E) containing C satisfies E.

Vertex Covers
In this section, we establish upper and lower bounds on (certain) vertex covers in a maximum degree d graph, which will be essential in proving the main result.
Proposition 2. Let C be a vertex cover of a connected graph G of maximum degree d. Then, |V | d+1 ≤ |C|. Proof. Let C be a vertex cover of a connected graph G of maximum degree d. Then, V \ C is an independent set. Since G is connected, each vertex in V \ C is incident to at least one edge with a vertex in C. Hence, there are at least |V \ C| edges between C and V \ C. Since each vertex in C has degree at most d, and we are done.
We now prove a lemma that is a key technical ingredient for the proof of our main result.
Theorem 2. Let d ≥ 3, let G = (V, E) be a graph of maximum degree d, and let I be a family of subsets of V . If I covers the vertex covers of G, then The above statement has been already proved by Razgon [Raz14b, Theorem 4], but here we give a new and shorter proof.
Our proof relies on the following theorem. Let G = (V, E) be a graph and let S ⊆ V . We denote by Theorem 3. Let G = (V, E) be a graph, S ⊆ V , and s ∈ S. Then, We define a function f s : If e contains s, then e is covered by a vertex in N s ⊆ D. Otherwise, if e does not contain s, then e is covered by a vertex in C \ {s} ⊆ D.
We now prove that f s is injective, so that from which the statement follows.
Proof. Let G = (V, E) be a graph of degree d, and let S ⊆ V . By induction on |S| ≥ 0, we prove that The base case S = ∅ is trivial. Let S = {s 1 , . . . , s n , s n+1 } where n ≥ 0. By the induction hypothesis, and we are done.
If t = min{|S| : S ∈ I}, then by Theorem 1 we have that, for all S ∈ I, Using the previous inequality, we get and we are done.
In establishing [Raz14b, Theorem 4], Razgon essentially proves that , which is slightly improved by our Theorem 3. As a consequence, also the bound in Theorem 2 slightly improves the 1 d 2 +1 log 2 1 + 1 2 d+1 −1 dependence on the maximum degree d in [Raz14b,Theorem 4]. Finally we observe that the bound in Theorem 3 is tight. It is easy to see that the vertex covers of the star of the degree d and center c, that is the graph whose edge set is {{c, r 1 }, . . . , {c, r d }}, are exactly {r 1 , . . . , r d } and {c} ∪ S for all S ⊆ {r 1 , . . . , r d }; hence R d has exactly 1 + 2 d vertex covers and exactly 2 d vertex covers containing c. If the graph G is the disjoint union of some stars of degree d, and S contains the center of each star, then G has exactly (1 + 2 d ) |S| vertex covers and exactly 2 d|S| vertex covers containing S.

Nice DNNFs
In this section, we observe that a DNNF computing a graph CNFs is equivalent to a DNNF that is nice, in that its internal gates have fanin at most 2, and its input gates are labelled by positive literals (Corollary 2). We also observe that nice DNNFs have nice certificates, in the sense that their underlying graphs are rooted binary trees whose leaves are labelled by positive literals (Proposition 6).
Proposition 3. Let D be a DNNF. There exists a DNNF D ′ equivalent to D, having fanin 2, and such that size(D ′ ) ≤ 2 · size(D).
Proof. Let D be a DNNF. An NNF D ′ equivalent to D and having fanin 2 is obtained by editing D as follows, until no gate of fanin larger than 2 exists: Let v be a ∧-gate with input wires from gates v 1 , . . . , v i with i > 2; delete the wires (v j , v) for j ∈ {2, . . . , i}; create a fresh ∧-gate w, and the wires (w, v) and (v j , w) for j ∈ {2, . . . , i}. The case where v is an ∨-gate is similar.
It is readily verified that D ′ is decomposable. Moreover, each wire in D is processed at most once (when it is an input wire of a gate having fanin larger than 2) and it generates at most 2 wires in D ′ , hence the size of D ′ is at most twice the size of D.
. An NNF is called negation free if its labels do not contain negative literals. Proof. Suppose D contains a gate u labelled with literal ¬x. Let D ′ be the DNNF on vars(D) obtained from D by relabelling u with the constant 1. We claim that an assignment satisfies D if and only if it satisfies D ′ . Let f be a satisfying assignment of D. By Proposition 1 there is a certificate T of D such that f satisfies T . We obtain a certificate T ′ of D ′ by relabelling u with the constant 1 (if u appears in T ). It is straightforward to verify that f is a satisfying assignment of T ′ . We apply Proposition 1 once more to conclude that f must be a satisfying assignment of D ′ . For the converse, let f be a satisfying assignment of D ′ , and let T ′ be a certificate of D ′ such that f satisfies T ′ . If T ′ does not contain the gate u then T ′ is also a certificate of D, and f is a satisfying assignment of D by Proposition 1. Otherwise, we obtain a certificate T of D from T ′ by relabelling u with the literal ¬x. Let f ′ be the assignment such that f (y) = f ′ (y) for all y ∈ vars(D) \ {x}, and such that f ′ (x) = 0. Since D is decomposable and T contains the gate u labelled with ¬x, no node of T can be labelled with the literal x. Thus T ′ cannot contain such a gate either, and f ′ satisfies T ′ . Since ¬x evaluates to 1 under f ′ , the certificate T is satisfied by f ′ as well. By Proposition 1, the assignment f ′ is a satisfying assignment of D. Because the function computed by D is monotone, we conclude that f must satisfy D as well. Clearly D and D ′ have the same size and maximum fanin. It follows that the desired negation free DNNF can be obtained by replacing every negative literal by the constant 1 in the labels of D. Corollary 2. Let E be a graph CNF and let D be a DNNF equivalent to E. There exists a fanin 2, constant free, and negation free DNNF D ′ equivalent to D such that size(D ′ ) ≤ 2 · size(D).
Proof. Let E be a graph CNF and let D be a DNNF equivalent to E. By Proposition 3, there exists a DNNF D 1 , equivalent to D, having fanin 2 whose size is at most twice the size of D. Since E is a monotone Boolean function, by Proposition 4, there exists a fanin 2 and negation free DNNF D 2 , equivalent to D 1 , whose size is at most the size of D 1 . Since E is a non constant Boolean function, by Proposition 5, there exists a fanin 2, negation free, and constant free DNNF D 3 , equivalent to D 2 , whose size is at most the size of D 3 . Let D ′ = D 3 . Then, D ′ is a fanin 2, constant free, and negation free DNNF, equivalent to D, whose size is at most twice the size of D.
Proposition 6. Let D be a DNNF and let T ∈ cert(D).
• If D is constant free, then the underlying graph of T is a tree.
• If D is constant free and fanin 2, then the tree underlying T has maximum degree 3.
• If D is constant free and negation free, then atoms(T ) ⊆ X.
Proof. We prove the first item (the rest is clear). Assume that the underlying graph of T is cyclic, so that in the underlying DAG there exist two distinct nodes v and w in T and two arc disjoint directed paths from v to w; in particular, w has at least two ingoing arcs in T , hence by construction w is a ∧-gate in D. By decomposability, no variables occur as labels of input gates in sub(D, v), which is impossible since D is constant free.

DNNFs Models
In this section, we prove two lemmas to which we critically appeal in the proof of the main result.
The first lemma (Lemma 1) allows to derive, given a DNNF D and a gate v in D, a DNNF D ′ whose models are exactly the models of D minus those models that correspond (in the sense of Proposition 1) to certificates of D using the gate v.
Let D be a constant free DNNF and let v ∈ D \ {output(D)}. Let elimgate(D, v) denote the DNNF obtained from D as follows. Let paths(D, v) be the set of directed paths from v to output(D) in the DAG underlying D; display each P in paths(D, v) by (w 0 , w 1 , . . . , w p ), where w 0 = v, w p = output(D), and p ≥ 1.
If there exists P = (w 0 , w 1 , . . . , w p ) in paths(D, v) such that w j is an ∧-gate for every j ∈ [p] = {1, . . . , p}, then elimgate(D, v) is the DNNF formed by a single 0-gate (computing the constant 0). Otherwise, for every P = (w 0 , w 1 , . . . , w p ) in paths(D, v) there exists j ∈ [p] such that w j is an ∨-gate; let e P = (w j−1 , w j ) denote the wire in P such that j is the smallest index in [p] where w j is an ∨-gate. Then, elimgate(D, v) is obtained from D by deleting the wires {e P : P ∈ paths(D, v)} in D. The process might create 0-gates in the form of fanin-0 ∨-gates, and undesignated sinks. Proof. Let f be a satisfying assignment of D such that f does not satisfy any certificate T of D with v ∈ T . Let T be a certificate of D such that f satisfies T ; such a certificate exists by Proposition 1. Clearly, v ∈ T . We claim that all wires in T are kept in passing from D to elimgate(D, v); therefore, T is a certificate of elimgate (D, v), and since f satisfies T , by Proposition 1, it also holds that f satisfies elimgate(D, v).
We prove the claim. Assume for a contradiction that a wire (u ′ , u) of T is deleted in passing from D to elimgate(D, v). Then, by definition of elimgate(D, v), there is a path from v to u ′ whose nodes are all ∧-gates. Since (u ′ , u) is a wire of T and thus u ′ is a gate of T , we have, by definition of certificates, that v is in T , which is a contradiction.
Conversely, let f be a satisfying assignment of elimgate(D, v). By Proposition 1, there exists a certificate T of elimgate(D, v) such that f satisfies T ; then, T does not contain 0-gates. In particular, T does not contain 0-gates that arise, in passing from D to elimgate(D, v), as fanin-0 ∨-gates on some directed path in D from v to output(D). Therefore, every gate w ∈ T is either a ∧-gate in D, so that the gate itself and all its input wires are kept in elimgate(D, v) by construction, or is a ∨-gate in D for which at least one of its input wires is kept in elimgate(D, v), and is used by T . Hence, T is a certificate for D. The node v does not have a directed path to output(D) in elimgate(D, v), because at least one wire is deleted on each such path in passing from D to elimgate(D, v); hence, v ∈ T . Since f satisfies T , by Proposition 1, f satisfies D.
The second lemma (Lemma 2) essentially establishes that if the models of a DNNF D computing a graph CNF have to satisfy a clause {x, x ′ }, and for some gate v in D it is the case that an input gate labelled x has a directed path to v in the DAG underlying D but no input gate labelled x ′ has such a directed path to v, then all certificates of D using the gate v hit x or all certificates of D using the gate v hit x ′ .
First, we observe that x ∈ atoms(T ) \ atoms(sub(T, v)). Otherwise, assume that x ∈ atoms(T ) \ atoms(sub(T, v)), and let w be the gate of D lying in inputs(T ) \ inputs(sub(T, v)) labelled x. By hypothesis, x ∈ atoms(sub(D, v)). Let w ′ be the gate of D in inputs(sub (D, v)) labelled x. Then, it is possible to construct a certificate S for the constant free DNNF sub(D, v) such that the gate w ′ is an input gate of S, and therefore x ∈ atoms(S); the underlying graph of S is a tree by Proposition 6. Now, the tree R obtained by replacing in T the subtree rooted at v by S is a certificate of D; in R, there exists an internal gate u such that both w and w ′ have a directed path to u in R; hence, u has fanin at least 2, so that u is a ∧-gate in D, which contradicts the decomposability of D (at u by variable x).
Second, since atoms(sub(T ′ , v)) ⊆ atoms(sub(D, v)), and x ′ ∈ atoms(sub(D, v)) by hypothesis, it holds that x ′ ∈ atoms(sub(T ′ , v)). Therefore, Now, the tree S obtained by replacing in T the subtree rooted at v by the subtree rooted at v in T ′ is a certificate of D; moreover, {x, x ′ } ∩ atoms(S) = ∅, contradicting the hypothesis that all certificates of D have a nonempty intersection with {x, x ′ }.

Lower Bound
In this section we prove the main result, an exponential separation between CNFs and DNNFs (Theorem 5).
Outline. We provide an informal outline of the proof. We separate CNFs and DNNFs by a class of graph CNFs naturally corresponding to a family of expander graphs (Theorem 1). Without loss of generality, the DNNF size of a graph CNF can be measured on DNNFs having fanin 2 and not containing negations or constants (Corollary 2); in particular, input gates are variables (or equivalently vertices of the underlying graphs). Let E be a graph CNF whose underlying graph (V, E) is an expander, and let D be a DNNF equivalent to E. In view of an application of Theorem 2, we construct a family of subsets of V covering the set of satisfying assignments of E, as follows. In each step, we pick a gate in D, and then we edit D so as to remove all certificates through that gate. As long as D has certificates, it follows from Proposition 2 that it must have at least |V |/(d + 1) distinct variables (vertices) labelling its input gates. We find a gate v in D such that the variables occurring below v in D are at least |V |/(d + 1) and at most |V |/2 (Lemma 3, Claim 1); this is sufficient and necessary to appeal to the expansion of the graph underlying E, and find in such a graph a matching of size proportional in |V | between variables occurring below v in D and variables not occurring below v in D (Lemma 3, Claim 2). We then argue using Lemma 2 that there is an equally large subset of vertices covering the satisfying assignments corresponding (in the sense of Proposition 1) to certificates through v (Lemma 3, Claim 3). We now remove the certificates through v from D, by following every path from v to the output gate of D until we find a wire that leads to an ∨-gate, and by deleting this wire. By Lemma 1, performing these deletions eliminates a satisfying assignment of D if and only if it corresponds to a certificate through v. The iteration terminates once the resulting DNNF is unsatisfiable. At termination, we have found a family of substes of V , each of size proportional in |V |, that covers all satisfying assignments of E; by Theorem 2, this family has size exponential in |V |. Since there is a bijection between the substes of V in the family and a subset of gates in the original DNNF D, the latter must have exponential size in |V | as well.
Proof. We now prove the main result. The following is the main technical lemma.
, and let D be a fanin 2, constant free, and negation free DNNF computing the CNF E. There exist K ⊆ D and a family I of subsets of vertices in G such that: (i) |K| ≥ |I|; (ii) I covers the vertex covers of G; (iii) min{|I| : I ∈ I} ≥ min{1,c} 2d 2 |V |.
Proof. We construct K ⊆ D and a family I of subsets of vertices in G as follows. Let D 0 = D, K 0 = ∅, and I 0 = ∅. For i ≥ 0, we define D i , K i , and I i maintaining the following invariants: (I1) D i is a fanin 2, constant free, and negation free nontrivial DNNF; At i = 0, the invariant (I1) holds by the hypotheses on D and the fact that the CNF E is neither unsatisfiable nor valid; the invariant (I2) holds trivially because K 0 = ∅. Let i ≥ 1, and assume that the invariant holds at i − 1.
Proof of Claim 1. Let (w 1 , . . . , w s ) be a directed path in D from an input gate w s to the output gate w 1 = output(D i−1 ) such that for all j ∈ [s − 1] if w j has input wires from gates w j+1 and w in D i−1 (by hypothesis, w j has fanin at most 2), then |atoms(sub(D i−1 , w))| ≤ |atoms(sub(D i−1 , w j ))|.
By (I1), it holds that D i−1 is nontrivial, so that there exists a satisfying assignment f ∈ sat(D i−1 ). By (I2), sat(D i−1 ) ⊆ sat(D). By Proposition 1, there exists a certificate T ∈ cert(D i−1 ) such that f extends T . Since D i−1 is constant free and negation free, atoms(T ) ⊆ V .
We observe that atoms(T ) contains a vertex cover of G. To see this, assume by way of contradiction that no vertex cover of G is contained in atoms(T ). Let g be the assignment defined by g(x) = 1 if and only if x ∈ atoms(T ). Since g extends T , we have that g ∈ sat(D i−1 ) by Proposition 1. Consequently, g ∈ sat(D). Since D computes E, it holds that g ∈ sat(D) if and only if {x : g(x) = 1} is a vertex cover of G. But by definition {x : g(x) = 1} ⊆ atoms(T ), a contradiction.
As G is a (c, d)-expander, G is connected and has maximum degree at most d, we have by Proposition 2, that any vertex cover of G has size at least |V | d+1 . Hence, |atoms(T )| ≥ |V | d+1 . Therefore, we have that If we have |atoms(D i−1 )| ≤ |V | 2 , we may choose v i = w 1 = output(D i−1 ) and are done. Otherwise, where the first inequality holds as d ≥ 3. Hence, there exists an index j in [s − 1] such that by construction, where the last inequality holds as d ≥ 3. We choose v i = w j .
Let v i be the gate in D i−1 given by Claim 1. If elimconst(elimgate(D i−1 , v i )), that is, the result of first eliminating gate v i and then constants in D i−1 , computes the constant 0 function, then the construction terminates with K = K i−1 and I = I i−1 . Otherwise, put Note that D i satisfies (I1) by inspection of the elimination operation. Moreover, D i satisfies (I2) as follows: Claim 2. There exists a matching M i in G between atoms(sub(D i , v i )) and V \ atoms(sub(D i , v i )) such that Proof of Claim 2. Let S = atoms(sub(D i , v i )). Since |S| ≤ |V |/2 and G is a (c, d)-expander, by (1) we have that Recall that |S| ≥ |V |/(d + 1) by Claim 1. We construct a matching M i between between S and V \ S in G as follows. Pick an edge {v, w} ∈ E with v ∈ S and w ∈ N S ⊆ V \ S; add {v, w} to M i ; delete v from S, w from N S , the vertices in S with no neighbors in N S after the deletion of w, and the vertices in N S with no neighbors in S after the deletion of v; iterate on the updated S and N S , until either S = ∅ or N S = ∅. At each step, we delete at most d vertices in S and at most d vertices in N S .
Hence, we iterate for at least 2d 2 |V |, and we are done.
Let M i be the matching in G given by Claim 2.
Claim 3. There exists a subset I i of vertices in G such that |I i | = |M i | and, for all certificates T of D such that v i ∈ T , and all vertex covers C of G such that atoms(T ) ⊆ C, it holds that I i ⊆ C.
Proof of Claim 3. We construct a set I i ⊆ V as follows. Let {x, x ′ } ∈ M i . By construction, x ∈ S = atoms(sub(D i , v i )) and x ′ ∈ N S . Since M i ⊆ E and D computes E, it holds that {x, x ′ } ∩ atoms(T ) = ∅ for all T ∈ cert(D); otherwise, by Proposition 1, an assignment extending T and sending both x and x ′ to 0 would satisfy E, contradiction. By the invariant, any certificate of D i is also a certificate of D; hence {x, x ′ } ∩ atoms(T ) = ∅ for all T ∈ cert(D i ). Then, by applying Lemma 2 to D i , x ∈ atoms(T ) for all T ∈ cert(D i ) such that v i ∈ T , or x ′ ∈ atoms(T ) for all T ∈ cert(D i ) such that v i ∈ T . If x ∈ atoms(T ) for all T ∈ cert(D i ) such that v i ∈ T , then add x to I i ; otherwise, add x ′ to I i .
By construction, I i contains exactly one vertex for each edge in M i , hence |I i | = |M i |. Finally, let T be a certificate for D such that v i ∈ T . We observe that I i ⊆ atoms(T ). Indeed, let w ∈ I i , and let {w, w ′ } ∈ M i be the corresponding edge in M i . By construction, all certificates containing v i have w as an atom; in particular, w ∈ atoms(T ). We conclude that, for all vertex covers C of G such that atoms(T ) ⊆ C, it holds that I i ⊆ C. Put: We conclude checking that the sets K = {v 1 , . . . , v L } and I, returned by the above construction, satisfy (i)-(iii) as required.
(i) The mapping sending v i ∈ K to I i ∈ I (i ∈ [L]) is, by construction, a surjective function from K to I. Hence, |K| ≥ |I|.
(ii) Let C be any vertex cover of G. As D computes E, the assignment f defined by f (v) = 1 if and only if v ∈ C satisfies D, that is, f ∈ sat(D). By (I2), there exists i ∈ [L] such that f ∈ sat(D i−1 ) and f ∈ sat(D i ). By construction, there exists a certificate T ∈ cert(D) such that v i ∈ T and f ∈ sat(T ); that is, atoms(T ) ⊆ {v ∈ V : f (v) = 1} = C. So, by Claim 3, I i ⊆ C. Thus every vertex cover of G contains some set in I; that is, I covers the vertex covers of G.  Proof. First note that E = ∅, because |V | ≥ 2 and (c, d)-expanders are connected. Let D be a DNNF equivalent to the graph CNF E. By Corollary 2, there exists a fanin 2, constant free, and negation free DNNF D ′ equivalent to D such that 2·size(D) ≥ size(D ′ ). By Lemma 3, there exist K ⊆ D ′ and a family I of subsets of vertices in G such that |K| ≥ |I|, I covers the vertex covers of G, and min{|I| : I ∈ I} ≥ min{1, c}|V |/2d 2 . Putting h(c, d) = min{1,c}·f (d) Now observe that |E| ≤ d|V |, because G has maximum degree d; thus the CNF E has at most d|V | clauses, each of at most 2 literals, so that d|V | + 2d|V | = 3d|V | ≥ size(E). Therefore, where g(c, d) is as in the statement, and we conclude that DNNF(E) ≥ 2 g(c,d)·size(E)−1 .
Theorem 5. There exist a class C of CNF formulas and a constant c > 0 such that DNNF(F ) ≥ 2 c·size(F ) for each formula F ∈ C. Indeed, C is a class of read 3 times monotone 2-CNFs.
Proof. By Theorem 1, there exists a family G = {G i = (V i , E i ) : i ∈ N} of (e, 3)-expander graphs such that |V i | ≥ 2 for all i ∈ N and |V i | → ∞ as i → ∞ (e > 0). Every graph in G is connected; in particular, it does not contain isolated vertices. Therefore E i is a CNF for every i ∈ N, and indeed E i is a read 3 times monotone 2-CNF.
Since |V i | → ∞ as i → ∞ and each graph in G satisfies (1), there exists an infinite subset I ⊆ N such that size(E i ) < size(E i+1 ) for all i ∈ I. Choose c > 0 and j ∈ I large enough such that g(e, 3) · size(E j ) − 1 ≥ c · size(E j ), where g(·, ·) is as in the statement of Theorem 4. It follows from Theorem 4 that DNNF(E j ) ≥ 2 g(e,3)·size(Ej)−1 ≥ 2 c·size(Ej) ; we take C = {E i : i ∈ I and i ≥ j}, and the statement is proved.

Corollaries
In this section we will prove the corollaries of Theorem 5 we sketched in the introduction.
Let F be a CNF-formula. We say that a clause C is entailed by F if every satisfying assignment of F also satisfies C. We say that a clause C ′ subsumes C, if C ′ ⊆ C. A CNF F is in prime implicates form (short PI) if every clause that is entailed by F is subsumed by a clause that appears in F and no clause in F is subsumed by another. Note that CNFs in PI form can express all Boolean functions but it is known that encoding in PI form may generally be exponentially bigger than general CNF [DM02].
Lemma 4. Every monotone 2-CNF formula is in PI form.
Proof. Let F be a monotone 2-CNF formula. We first note that trivially no clause in a 2-CNF subsumes another. Now let C be a clause entailed by F and assume by way of contradiction that C is not subsumed by any clause of F , that is, every clause in F contains a positive literal not in C. Let C ′ be the clause we get from C by deleting all negative literals. We claim that C ′ is entailed by F . To see this, consider a satisfying assignment f of F . Let f ′ be the assignment we get from f by setting the variables that are negated in C to 1. Since F is monotone, this is still a satisfying assignment of F and thus of C. Consequently, C is satisfied by one of its positive literals in f ′ and thus in f . Thus f satisfies C ′ and it follows that F entails C ′ . Now let f be the assignment that sets all variables in C ′ to 0 and all other variables to 1. Since C and thus also C ′ is not subsumed by any clause of F , the assignment f satisfies F . But by construction f does not satisfy C ′ which is a contradiction.
Remember that the formulas of Theorem 5 are monotone 2-CNF formulas. Then we directly get the promised separation from Lemma 4 and Theorem 5.
Corollary 3. There exist a class C of CNFs in PI form and a constant c > 0 such that, for every formula F in C, every DNNF equivalent to F has size at least 2 c·size(F ) .
It follows that L 1 can be exponentially more succinct than L 2 for any two representation languages L 1 ⊇ PI and L 2 ⊆ DNNF. In particular, this holds if L 1 ∈ {PI, CNF, NNF} and L 2 ∈ {d-DNNF, DNNF}. Here, d-DNNF denotes the language of deterministic DNNFs, that is, DNNFs where subcircuits leading into DNNFs where subcircuits leading into a common ∨-gate never simulataneously evaluate to 1. This answers several questions concerning the relative succinctness of common representation languages [DM02].
We also observe that DNNFs are not closed under negation.
Lemma 5. There exist a class D of 2-DNF formulas and a constant c > 0 such that, for every formula D in D, every DNNF equivalent to ¬D has size at least 2 c·size(D) .
Proof. Let C be the class of 2-CNFs from Theorem 5. Let D be the class of 2-DNFs we get by negating the formulas in C. Now negating D gives the class C again, for which we have the lower bound from Theorem 5.
Observing that DNF is a restricted form of DNNF, we get the following non-closure result which was only known conditionally before.
Corollary 4. There exist a class D of DNNFs and a constant c > 0 such that, for every formula D in D, every DNNF equivalent to ¬D has size at least 2 c·size(D) .

Conclusion
We proved an unconditional exponential separation between the representational power of CNFs and that of DNNFs and discussed its consequences in the area of knowledge compilation [DM02]. Let us close by mentioning directions for future research.
In order to prove the lower bound of Theorem 5, we generalized arguments concerning paths in branching programs to the tree-shaped certificates of DNNFs. It would be interesting to know whether other lower bounds for branching programs can be lifted to (suitably restricted) versions of DNNFs along similar lines.
Recent progress notwithstanding [BLRS13], several separations between well-known representation languages are known to hold only conditionally [DM02]. For instance, it is known that DNFs cannot be compiled efficiently into so-called deterministic DNNFs unless the polynomial hierarchy collapses [SK96,CDLS02]. It would be interesting to show this separation unconditionally.
Finally, there is a long line of research proving upper bounds for DNNFs and restrictions (see [Raz14a,OD14a,OD14b] for some recent contributions). We believe that these results should be complemented by lower bounds as in [Raz14a,Raz14b], and hope that the ideas developed in this paper will contribute to this project.


A hierarchical map building for SLAM used in ruins environments

As to the morphological characteristics of the seismic damage of the interior ruins, this paper presents a map building mechanism by using the analytic hierarchy process and proposes a hierarchical simultaneous localization and mapping (SLAM) algorithm that is based on the hybrid metric-topological map. The architecture of the ruins-oriented SLAM system consists of three layers. In the intermediate layer, the global environment is partitioned on the basis of the spectral clustering that can reduce the computational cost and ensure the global consistency. A semantic map is built in the top layer where the node identification and logic location are included. According to the characteristics of the seismic destruction, an aggregation degree of features is proposed as an identification method for building topological nodes to improve the interactivity. In the bottom layer, details of the morphology of the destruction are described as grid maps that ensure the environmental adaptability. Through experiments, the ability of environmental description and the availability at an artificial ruins environment are verified.

I. INTRODUCTION
The frequent outbreak of natural disasters, such as earthquake, is a serious problem in recent years [1]. The development of autonomous robots for the search and rescue then gradually becomes a hot topic. For a rescue robot, the ability of the localization and mapping is a precondition. The internal environment of ruins which is formed by random seismic damage is complex and unknown. The process of localization and mapping lacks an absolute reference. That is both is highly relevant and constantly iteration, which is the problem of simultaneous localization and mapping (SLAM). SLAM [2] addresses a problem how to acquire the localization of a mobile robot while simultaneously building a spatial map of environment without information about the global position (e.g. GPS).
The research of SLAM gradually extends from the theoretical simulation to the practical application. In terms of the search and rescue after earthquakes, Behdad Soleimani [3] at Tehran University solved the SLAM problems at the external of ruins via the extraction of invariant features from aerial images. To realize the localization and mapping of the rough terrain in narrow canals, a three-dimensional map representation called S-DEM [4] was proposed by Keiji Nagatani. In [5], Alexander Kleiner presented the DSLAM algorithm based on the method of non-selfish sharing of information to execute the rapid rescue while the communication is intermittent.
The above methods, however, did not consider the environmental model based on morphological characteristics of the seismic damage of the interior ruins. As a model of the environment, the map building not only desires an appropriate description for the target environment, but also takes into account the computational complexity of filters.
In this paper, by using the analytic hierarchy process (AHP), a hierarchical SLAM based on the hybrid metrictopological map is proposed. To begin with the division of the environment based on the spectral clustering method in the intermediate layer, the computational efficiency and global consistency are balanced. Then on this basis, details of seismic damage are described as grid maps in the metric layer. In the semantic layer, the damaged structure is recognized by the aggregation degree of corner features to build a semantic map. Finally, the simultaneous localization and mapping of the internal ruins is achieved at different hierarchies.
The rest of this paper is organized as follows. Section II presents the modeling of the internal ruins environments. Section III describes the architecture of the hierarchical SLAM that is based on the hybrid metric-topological map. The implementation of the metric map on bottom layer and the semantic map on the top layer are presented in section IV. Section V presents system experiments and analysis. Finally, some conclusions are summarized.

II. DECISION MODEL OF THE MAP REPRESENTATION OF INTERNAL RUINS ENVIRONMENTS
As one of core issues of SLAM, the map representation determines the estimation method and the data association. A map model for the ruins environment needs to describe the seismic damage inside buildings and simultaneously meets needs of the rescue mission. Based on the analytic hierarchy process (AHP), a building mechanism of the map 978-1-4799-5825-2/14/$31.00 ©2014 IEEE Proceeding of the 11th World Congress on Intelligent Control and Automation Shenyang, China, June 29 -July 4 2014 representation for the interior ruins is analyzed. Then on this basis, an optimal decision of the environment model can be obtained.

A. Analysis of the Irregular Seismic Damage inside Semi-Collapsed Buildings
In this paper, the ruins environment means seismic disasters. The internal environment is composed of severely damaged or completely destroyed building components [6]. The seismic damages are related to the structural type and construction time of buildings, and the intensity of earthquakes. Degrees of seismic damages [7] are divided into five grades. The grades of basically intact, minor damage, moderate damage and seriously destroy are mainly considered in this paper. Among them, the damage degree of bearing members is an especially important criterion to the grade evaluation of the seismic damage of buildings. Therefore, the environmental modeling needs to focus on the seismic damage of the columns and the construction walls. The schematic diagrams of the four destructive grades of building components are shown in Table I. Damaged structures, such as cracks in walls and columns or collapsed deposits, partially appear inside the semi-collapse building, as shown in Fig. 1. The randomized and uneven distribution of damaged structures appears irregular boundary on a two-dimensional plane. Such a complex environmental characteristic requires an appropriate map representation for the SLAM problem.
Besides, because of the large-scale search and urgent time, a SLAM algorithm for the auxiliary rescue or the disaster evaluation should has strong computational capabilities to achieve the rapid localization and mapping.

B. The Decision Model of Map Representation Based on Analytic Hierarchy Process
On the basis of AHP algorithm [8], a building mechanism for the decision model of map representation needs to determine relative weights of the decision criteria and relative rankings of alternatives. The building process is shown in  The model criterion layer is created in accordance with created targets and evaluation criteria of SLAM problem. This layer includes applicability C 1 , accuracy C 2 , computational cost C 3 , and interactivity C 4 . Map representations at the scheme layer are topological map S 1 , metric map S 2 and hierarchical map S 3 . Please refer to [9] for a detailed description of the comparison of alternative map representations.
2) Calculation of Judgment Matrixes and the Single Hierarchical Ranking: The influence from each item of criterion layer to the target is calculated and pairwise compared. And then, the proportion of each criterion is determined to sort the importance.
According the aforementioned analysis of the seismic damages, the pairwise comparisons are calculated to construct a judgment matrix A=(a ij ) n×n . The maximum eigenvalue and eigenvector of A are calculated. And then, the consistency ratio (C.R.) is checked based on C.R.=C.I./R.I.=0.0564<0.1. Among them, C.I. is the consistency index and R.I. is the random index.
Thus, the resulting eigenvector can be the weight vector which means the single hierarchical ranking, as shown in Table II. ω C ={ω C1 , ω C2 , ω C3 , ω C4 }={0.4369, 0.0467, 0.4066, 0.1099}. The single hierarchical ranking of the scheme layer and the criterion layer is calculated in the same way. The weight C S ω is shown in Table III.

3) Decision of the environment model in ruins:
According to the resulting ranking of criterion layer and scheme layer, the weight ω S between each scheme to the target can be calculated based on (1).
The weight of scheme S 3 is the highest. Thus, a hierarchical map representation is the optimal decision of the environmental model. According to the above analysis, a theoretical basis is provided for the hierarchical map representation of SLAM in ruins environments.

III. THE OVERVIEW STRUCTURE OF THE HIERARCHICAL SLAM IN RUINS ENVIRONMENTS
A. The Hierarchy of the SLAM System According to the above decision model of the map representation, the architecture of the ruins-oriented SLAM system based on a hierarchical map is shown in Fig. 4. The architecture is composed of three parts: semantic layer, intermediate layer and metric layer. The intermediate layer includes several sub maps which are obtained by dividing a large scale environment. On each sub map, the damaged structure is identified as a topological node which represents the distribution of damage of walls or columns in the semantic layer. Environmental details inside each sub map are described by occupancy grids.  The building flow of the three layers is shown in Fig. 5. The process is started with the preprocessing of original data which is relative distances from a laser range finder and control instructions of a mobile robot in this paper. After the extraction and estimation, positions of corner features and the robot are achieved. Then, the current environment can be divided into several sub maps based on the spectral clustering. The detail of the sub map is represented as grids by updating the occupy probability in the underlying metric layer. And in the upper semantic layer, on the basis of calculation of the aggregation degree, the damaged structure is identified at each sub map.

B. Environmental partitioning based on spectral clustering
Based on the map segmentation, a large-scale environment is divided into several local regions with an appropriate scale. Through the segmentation, each local region is independent, and therefore the computational cost is reduced. Unfortunately, the inherent connection between local maps is also cut off. Discarding the interregional correlation would affect the accuracy and consistency of SLAM. A pivotal issue is how to balance the computational complexity and the global consistency during the map segmentation.
An environmental partition based on spectral clustering [10][11] is utilized in this paper. According to the similarity of environmental structure, a recursive normalized cut is proposed to obtain local sub maps.
This approach considers each robot's pose as a node of an incrementally built graph whose edge represents the sensedspace overlap (SSO) between two observations which have sensed at each node. The SSO is a pairwise similarity measure which reflects the common space between observations, is defined as: Where M(z t ) stands for the set of map elements sensed in a given observation z t . The map segmentation is then abstracted as an observation graph cut issue.
Finding the exact min-Ncut bisection is a NP-complete problem, thus we use an approximate approach which analysis the spectrum of the Laplacian matrix of the graph based on spectral bisection of graphs. Please refer to [12] for a detail description of spectral bisection of a graph. The whole process of environmental partitioning is summarized as shown in Fig. 6.

IV. IMPLEMENTATION METHODS OF THE HIERARCHICAL MAP BUILDING
A. Topological map building and localization in the top layer 1) Definition of topological nodes: Random earthquakes cause serious damage to the building. The existing methods of node identification are mainly used to identify specific locations of the structured indoor environment, such as walls and doors. Thus, those methods cannot be applied in ruins directly. For missions of rescue and assessment, the damage of bearing members is a main target area. Therefore, we define walls and columns with cracks or breakages and collapsed deposits as topological nodes to construct the semantic space.
Among them, the topological node N is refers to the seismic damage region R i which includes walls and columns with cracks W j and collapsed deposits D k .
2) Identification of topological nodes: In a twodimensional mapping plane, general destructive patterns appear complex boundaries which can be abstracted as a plurality of intersecting line segments. We identify the topological node by the aggregation degree of corner features (intersection point of lines).
Firstly, the iterative end point fit (IEPF) [13] is used to extract the line segment in polar coordinates. The process of feature extraction is expressed as: Wherein, Z 0 is the original observation. l i is the line segments feature and g j is the corner feature. S g is the set of extraction results of corners. Secondly, the aggregation degree of corner feature is used to identify the damaged structure. Density [14] of corner features in the detected area is calculated based on (6).
Wherein, (g) i r refers to the weight which is proportional to the distance from the corner. σis the integration degree of each point.
Special irregular contour of the damaged structure forms a region of dense corner features. Density f(p) of each corner feature is calculated. The function value is proportional to the aggregation degree of corners. The feature which has a greater value than the preset threshold can be joined into the set of damaged structures.

B. Metric map building and localization in the bottom layer
Grid maps which have no specific limitations of the environmental type are used to describe details of seismic damage morphology and ensure the accuracy.
Bayesian estimation is used to iteratively update the probability of occupancy grids. The initial state of grids is determined according to the expression range and detection accuracy of the laser range finder [15]. The state of grid is updated by "odds" on the updated phase. And the formula of iteratively update is expressed as: Among them, ( ) , | ,

V. EXPERIMENTS AND ANALYSIS
The effectiveness of the proposed hierarchical SLAM is verified in two environments. As shown in Fig. 7, a shape shifting search and rescue robot AMOEBA-I which carries a laser scanning distance sensor URG-04LX is used in these experiments. The sensor is carried on the middle module of the robot. The height of sensor from ground is about 35 cm. The distance range of the measure is 4 m and the angular range is 240 degree. The ability of the hierarchical map model to describe irregular obstacles is verified in an office environment. This is a cluttered indoor environment. Legs of chairs and tables, books and cardboard boxes on the ground are randomly distributed which appear complex contours on the observation plane. The result of SLAM algorithm on the bottom metric layer is shown in Fig. 8. According to partial enlarged views, grid maps can accurately depict the complicated shape that means it can completely describe irregular obstacles within a clutter environment. Therefore, the map representation is appropriate for the expression of the interior of ruins environment.

B. Experiment II
In order to simulate the afflicted scene after earthquakes, an artificial ruins experimental environment is constructed according to the seismic damage morphology. The correspondence between damage structures is shown in Fig. 9. Crevices in the wall are simulated by gaps between cartons. And the broken building components are randomly stacked on both sides of the corridor to simulate the local collapsed deposits. The experimental result illustrates that the algorithm can accurately describe the environmental structure. In this paper, a hierarchical map is used for the robot environmental modeling. The semantic map provides the distribution of damaged structures for the less complexity. A grid map is used in each partition to represent details of the irregular obstacles. On the premise of ensuring the computing power, a composite map is built to identify and describe the damaged structure.

VI. CONCLUSION AND FUTURE WORKS
Facing the search and rescue after earthquakes, a decision model of the map representation based on analytic hierarchy process has been proposed, according to the random seismic characteristics of the bearing members inside a semi-collapsed building. A collection of environmental expression was built by using the proposed hierarchical representation based on the hybrid metric-topological map. The environment was divided into sub maps on the basis of spectral clustering algorithm. Then, damaged structures were identified as semantic nodes and seismic damaged shapes were described as occupancy grids in each partition. Experimental results in a cluttered indoor environment and an artificial ruins environment were presented to verify the efficiency and applicability of the algorithm. To further improve the reliability and robustness of SLAM for actual search and rescue, we will pay more attention to the extensive application of multi perception information such as visual information and 3 dimensional point cloud data in the future.


extfieldsofstudy:
Computer Science
Physics

text:
FDCHQHP:A Fortran Package for Heavy Quarkonium HadroProduction

FDCHQHP is a Fortran package to calculate the transverse momentum (p_t) distribution of yield and polarization for heavy quarkonium hadroproduction at next-to-leading-order (NLO) based on non-relativistic QCD(NRQCD) framework. It contains the complete color-singlet and color-octet intermediate states in present theoretical level, and is available to calculate different polarization parameters in different frames. As the LHC running now and in the future, it supplies a very useful tool to obtain theoretical prediction on the heavy quarkonium hadroproduction.

I. INTRODUCTION
From the discovery of heavy quarkonium in the 1970s, both the theoretical and experimental study on production and decay of J/ψ and Υ, plays an important role in the research on the nonperturbative and perturbative aspects of QCD. The factorization framework of non-relativistic QCD (NRQCD) had been proposed to study the production and decay of heavy quarkonium [1]. By extracting the NRQCD long-distance matrix elements (LDMEs) from the matching between theoretical prediction and experimental data, the NRQCD calculation had given a good description on the transverse momentum distribution (p t ) of heavy quarkonium production at hadron colliders. But the polarization of heavy quarkonium hadroproduction had been an open question for many years.
In the past, many works has been done by applying the Fortran codes generated by our Feynman Diagram Calculation(FDC) package [2] which is written in REDUCE and RLisp and can generate Fortran codes automatically. In the recent two works [3,4], we have presented the complete next-to-leading-order(NLO) predictions to the p t distribution of yield and polarization for J/ψ and Υ prompt hadroproduction. It took us a lot of efforts to manipulate running of more than ten thousands the numerical calculation jobs and control the final numerical precision when more than hundred numerical results from sub-processes must be sum up for each polarization result at each p t point. Finally, based on the collection of all the Fortran codes generated by using FDC, the efforts and experience on job submission and numerical precision control are implemented in an independent package FDCHQHP, it can be easily run by inexperienced user to obtain any desirable results on heavy quarkonium hadroproduction at any given experimental condition. It is well known that only prompt experimental measurements on heavy quarkonium hadroproduction are available up to now. Therefore, to fit the experimental measurement, the theoretical prediction must be on the prompt production. FDCHQHP supplys the complete calculation on the p t distribution of yield and polarization for prompt hadroproduction, which includes direct production part and feed-down contribution form excited states production.
From another aspect, as the LHC running at the current energy, the p t distribution of yield are measured up to large p t region (e.g. the measurement of Υ by ATLAS are up to 70GeV [5]), and there are more polarization parameters measurement e.g. λ θ , λ θφ , λ φ in different reference frame for heavy quarkonium polarization, and also the frame-independent polarization parameter λ [6]. What's more, in the next years, the LHC running energy will upgrade to 14 TeV, and there will be more high-precision measurements to be published in the future. Although we have already finished a few theoretical prediction in comparison with the experimental measurements in past two work. There are still a lot of experimental measurements in future need the theoretical predictions which needed large scale calculation in supercomputer environment. There are already many request for theoretical prediction in different experimental condition, we are ready to make the package public so that it can be used to obtained theoretical prediction by other people.
This package includes 6 channels and 76 sub-processes as summarized in above  [8,9] are included in the package for convenience. Please note that, if the file makefile was changed, you should execute make clean first, and then execute make, otherwise, just make is OK.

III. THE DIRECTORY STRUCTURES AND PROGRAMS
A.

Structure of directories
All the Fortran codes are collected in the folder FDCHQHP. The folder's structure is like : • source/<channels>/<sub-processes>/fort/ • fdc2.0/f77/src/src64d ifort(src64d ifort p100, src64q ifort)/ , respectively. The <sub-processes> under each <channel> are the ones listed in Tab. I, they are named with 'g'or 'q', which denotes the partons in initial and final states except the heavy quarkonium. For example, gggg denotes g + g → J/ψ(Υ) + g + g, qggq denotes g + q(q) → J/ψ(Υ) + g + q(q), and so on. Without specification of quadruple precision, all the calculation are double precision calculation. Especially, in the channel named 3pj 8, the sub-process ggg loop32, denotes the one-loop part for g + g → J/ψ(Υ) + g(one-loop) with quadruple precision calculation Fortran source involved. A corresponding sub-process, ggg other denotes the remaining part aside from the loop part(i.e. ggg loop32 ) with double precision calculation.
We know that the physical processes are divided into several channels, and a channel are divided into a few sub-processes. Therefore, the final result in physics should be the summation of all of them. For any sub-processes, the Fortran source codes and the corresponding compiled executable program (i.e. int) are placed in the fort directory.

B. Programs
In the fort directory, there are several kinds of files, and the most important one for users is input.dat which present all the input to run the Fortran program. We describe them in the following.

Introduction to the file input.dat
For the direct production part, including 1 S 8 0 , 3 P 8 J , 3 S 1 1 , 3 S 8 1 , as well as for the feed-down contributions from χ c in the directories chic-3pj 8 and chic-3s1 8, the file input.dat are completely same, just as Tab.II.
As an example, when we calculated the p t distribution of polarization parameter(α) for J/ψ prompt production under a certain experimental measurement (here denotes the rapidity cut) at the LHC [3], we prepared two input.dat files which uses the format showing in Tab.II with me2 type set as 4 for longitude polarization and 1 for transverse polarization; Other parameters can be easily changed according to the detailed description in Tab.II.

Introduction to other files
Users can skip over this block directly, if they do not concern the details of other files.
• The files prefixed with amp, are the amplitude calculation part the ampN.f is the amplitude for the Nth tree-diagram, amplN.f is that for the Nth loop-diagram, amplNn2.f, amplNn3.f and amplNn4.f are the 1 ε IR , 1 In amps2.f file, it calls all the diagrams contribution and squares the matrix elements, there are different blocks for the tree-diagram part, the loop-diagram part, the remaining term from divergence canceling part. It also contains the different choices for polarization direction and frame.
• The files prefixed with incl contain the common blocks and the definition of data, used in the Fortran programs.
• makefile, list the rules, libraries with location to link and compiler for building the executable file int.
• parameter.f and parameter1.f are two files defining the physical parameters in calculation.
• func.f, generate a certain phase space point with a set of random number and calculate the cross section at the point.
• int.f, the main program, read all input and calls the Monte Carlo program for phase space integration.
• genppp.f, and coupling.f calculate all the Lorentz invariable and constants in the programs.
• bases.data and the files prefixed with plot, are generated automatically in calculation, which are usually not used.
• input.dat, lists the input parameters and provides an interface for users.
We make two tables in this paper to explain it as Tab. II.
• The files fresult.dat and convergence.dat, are the results we expect after programs running. The former is used in the data processing, and the latter is the precision analysis for the integration. They will be replaced by new ones if the executable int runs again. We will give some detailed introduction in the section OUTPUT FILES.
• Other files not listed.

A. Calculation
As we know, a complete prompt result should be the sum up of all the channels with fixed LDMEs. Users can also do the calculation for one sub-process independently. The following is the calculation introduction.

calculation for a sub-process
Once the programs compiled successfully, the executable file int will be found in the fort folder. Two calculation methods are available: It may take a few seconds to a few hours depending on the process if just one CPU is used. Usually, the loop diagrams will take more time. No matter which method is used, the input.dat and cteq6m.tbl(or other PDF used) files should be placed in the calculating folders correctly.

calculation for a channel or channels
Go into the source directory, execute ./run-all to do the calculation for all channel. It may take users several hours to a few days depending on the calculation precision and computer.
By modifying the variable $op in line 4 with any channel name in the file run-all, such as op = ′ 3s1 1/ * /f ort ′ , users can obtain the calculation result for the specified channel.
In the script run-all, we have processed the data collecting step. The results for different channels are listed in the corresponding folders, named as data *.

V. LARGE SCALE CALCULATION IN SUPERCOMPUTERS
To simplify the tedious manual work, we have developed a bash-shell script package WBIN to do the numerical calculation automatically on supercomputers. This package provides a very simply interface and some commands for users. It has been realized to automatically parallize the calculation, generate and submit jobs, analyze the accuracy of the results, fix the suitable number of CPU to each job, separately calculate the positive and negative part of a sub-process when needed, realize the test calculation(for 3 p t points) and the real calculation(for all p t points), check the state of each job, and finally collect all the data. It is customized for our package, but also could be adapted in any similar application which use BASES MC integral program.

A. Installation
The package can be downloaded together with the Fortran package, named as WBIN.taz.
There are two steps to install it: • execute tar -xzvf WBIN.taz ; • execute cd WBIN/ • add the absolute directory of WBIN to the default search path, e.g. add the line 'set path=($path $HOME/WBIN)' in ∼/.cshrc .

B. Usage
In the calculation of transverse momentum distribution of yield and polarization parameters of heavy quarkonium hadroproduction, all the results used the p t as the their x-coordinate and there are usually many p t points. The result at these p t points can be calculated independently. Instead of the sequential calculation mode for all the p t points at our personal computer environment, the calculation could be parallelled for all he p t points at the same time in supercomputer environment. Therefore, the parallelization scale will become very large dramatically. Finally, there are many jobs have to be submitted to finished in the parallel computing environment. And each job is for calculation on one polarization component, one sub-process at one p t point with one experimental condition.
We need to do a three p t points test running with low accuracy requirement for the result of each channel which is obtained by summing up all the sub-directories in this channel. In the test calculation, the number of sample point in the Monte Carlo integration is adjusted to obtained the required accuracy, as well as the number of CPU cores is adjusted to keep the calculation finish at two hours. With all the information obtained in the successful test calculation, we can arrange how to run program in each sub-process so that the accuracy requirement for the final result of each channel can be achieved. For each jobs, the number of CUP cores and number of sample point in the Monte Carlo integration are fixed by applying the information obtained in test calculation so that every jobs will be finished within two hours by themselves. Then the results for all p t points can be obtained in the real calculation stage.

preparation
The most important preparation is the method of submitting jobs, which depends on the job management system in supercomputers. In WBIN package, two files refer to this method: the file Submit and sub job. It needs users themselves to understand the method of submitting jobs in their supercomputer environment to modify it properly. Here we provide two examples which was applied by us at the Shenteng7000 and TH-1A in China.
• For Shenteng7000, we give the example as files Submit shenteng and sub job shenteng.
It adapts the LSF(Load Sharing Facility) HPC7.0 job management system at the Shenteng7000. There is no upper limit for the total number of CPU cores applied by users, but every job submitted is added in task queues together with other users' jobs.
• For TH-1A, we give the example as files Submit th-1a and sub job th-1a.
In TH-1A, it adapts the method of batch jobs. We can apply a certain number of CPU nodes within a upper limit, in these CPU nodes, we can arrange our small scale jobs.
Once the CPU nodes are allocated to a user, the user's jobs will monopolize the CPU nodes until all the jobs are over. It is called a non-preemptive scheduling system.

test calculation
For the usage of WBIN package, these steps will be done by manual operation or by shell script as the following shows.
• A: Make sure all of the sub-processes in folder FDCHQHP/source have been compiled successfully.
• B: Go into the directory source, execute echo $PWD > dir s.
• C: Make a directory any where which is specially used to do the calculation, e.g.
• D: Copy the files all-pwd, dir s, input.dat, cteq6m.tbl in source folder into the folder work.
• E: Go into the folder work, and use the template file input.dat to generate more input files.
In details, copy input.dat to other files named as inputX.dat with X can be any letters except '.', and then modify inputX.dat with the options you want according to Tab. II. The number of inputX.dat is not limited, such as inputlp.dat, inputtp.dat, inputdirectlp.dat ....
• F: If users want to calculate part of the total 76 sub-processes, the user can just keep the corresponding lines and delete others in file all-pwd.
• G: Execute init test job. The flow for this step: 0) generate the folders redo job, done job in work ; 1) generate the directories according to inputX.dat, all the X value will be the directory name in the folder work ; 2) generate the multiple sub-directories in directory X according to the file all-pwd, such as X/1s0 8/gggg/, and copy the inputX.dat file into X/1s0 8/gggg/input.dat; 3) generate the sub-directories in X/1s0 8/gggg/ according to p t points in the current input.dat and name them as tedir1, tedir2, tedir3..., for test calculation, and then generate input.dat in tedir1, tedir2, tedir3..., with just one p t point as the p t list at the end of input.dat; 4) link the file cteq6m.tbl in the work folder to tedir1, tedir2, tedir3...,;
• H: Submit jobs with executing Submit.
It will submit all the jobs listed in the file work/redo job/job list. Each line in the file denotes a job and looks like:  e) The step J: is also working for real calculation.

data collection
The test and real calculations make up the whole process for the calculation in supercomputers. Finally, execute final data collect to collect the data for all the p t points in all the sub-processes, and then give the final results for the 6 channels.
C. some detailed processing

control the calculation precision
Using the Monte Carlo integral program, the error estimate and the square root of integral points are in inverse ratio. Therefore, if we know the error for a certain number of sample points, it's easy to obtain the error for another number of sample points.
where N is the total number of sub-processes in a channel, c i is the central value for each sub-process, e i is the absolute error for each sub-process, and r i is the relative error for each sub-process. Therefore R is the relative error for a channel. Consuming the maximal relative error for a channel we desired is R 0 (the default value is 1%), when R > R 0 , we can easily obtain that in each sub-process, the relative error should be Using the new error r ′ i , the error r i and the number of sample points Ncall in the last calculation, we can get the new number of sample points Ncall ′ : Ncall ′ is the new number of sample points in the next calculation process. It may do several iterations until R < R 0 is achieved. This flow chart of these steps is shown in Fig. 1.
The calculation includes two steps. The first step is test calculation, in which we set a less precise error value r i = 1% for each sub-process, then do the calculation only at three p t points which are the minimal, intermediate, and the maximal number of p t points listed in input.dat. In this step, the relative error R maybe not meet the precision requirement for the result of a channel in the case of large number cancellation among different sub-processes in this channel. Then another iteration will work until R is small enough just as mentioned above. The second step is real calculation, in which all the p t points listed in input.dat will be calculated. Using the method of linear continuation for the number of sample points at the 3 p t points calculated in test step, we can obtain the numbers of sample points for all other p t points listed in input.dat. In this case, the calculation precision in the whole p t region could be guaranteed R 0 .

other details
• In some sub-processes, for the p t points in the region where the central value changes from positive to negative or in inverse, the precision will get worse in such region. In this case, we adapt to do the calculation for positive and negative part, respectively, and then sum up these two parts as the central value for the sub-process.
• In each sub-process, the number of sample points decides the run-time, so it is easy to estimate the number of CPU cores and the time needed in large number of sample points through the information obtained in small sample number case.
• When using the same number of CPU cores, run-time for different sub-processes are quite different. So we need to allocate the suitable number of CPU cores for each job to keep it can be finished in two hours.
• In WBIN package, we decide the priority to submit job by the number of CPU cores needed and then estimated run-time for the jobs.

a brief guide for the shell scripts
You can skip this subsection if you just want to run the package. It is presented here just for user who want to involved in modification of program in the package.
• Initialize the directories and the job list in test calculationinit test job, get test pt • Job Controllingjob control • Submit the jobssub job • Analyze the running state of jobanly job • Analyze the precisionanly precision In FDCHQHP package , the result for a channel "n" are summing up of all the subprocesses in this channel as: where i, j, k denotes the partons, σ a n is the result in the sub-process (corresponding to a) directories without the suffix '-m' in folder names, and σ b n is the result in the sub-process directories with '-m' suffix in folder names. The second part in the Eq.(5) only appears in the channels 3pj 8 and chic-3pj 1. It depends on the choice of µ Λ which usually take the value of m Q . Users can take different µ Λ values to study the dependence on NRQCD scheme, as we did in Ref. [4].
Considering the kinematics effect for the feed-down contribution from high excited states, could be performed after obtaining the result of a whole channel.

VII. OUTPUT FILES AND THE PRECISION ANALYSIS
As mentioned above, the files fresult.dat and convergence.dat are the output files needed.
The format offresult.dat is shown in Fig. 3, and the format of convergence.dat are shown in Fig.4 and Fig.5.
The following are some notes for convergence.dat: • In Fig.4, the first two lines are the beam energies for colliders; the 3rd line is the p t point, the 22nd shows the number of sample points (given and real used).
• The 38th and 41st lines show the maximal times for the grid optimization step and integration step, and the 39th and 42nd lines give the highest accuracy we expected in the two steps, usually they are set to be 0.01%.
• In Fig.5, the lines 45th to 49th lines give a statistical run-time for the calculation on one p t point. The lines from 60th to 69th are the result of each iteration in the grid optimization step. From left to right, they are the iteration number, the efficiency, the portion of the negative value, integration value, the relative error, the cumulative result for the data on the front(the absolute error), the total relative error, the time when each iteration finished.
• Lines 78th to 87th, are the reliable result in the integration step. In this step, after the grid optimization step, the results for each iteration are close to each other, and the final cumulative central value is the average of these data, the corresponding relative error(e.g. the 8 column of line 78, which is 0.077%) is the 1/ √ 10 times as the first line(e.g. the 8th column of line 87, 0.244% ).
There are two standards to ensure that the convergence behavior of integration is good enough: 1), the relative error in the last 10 lines of convergence.dat is average and small enough, no more than 1% (if the accuracy is as high as almost 0.01%, the less average result also could be accepted); 2), on the base of 1) achieved, the final cumulative accuracy in the last line, which is also the accuracy presented in fresult.dat in the form of absolute error, should meet the required precision when all the sub-processes are summed, since we know that large number cancellation usually appears in the summation of the sub-processes.
In fact, we also require that the relative error for each sub-process itself is no more than 1%( when it is larger than 1% , the relative error in each iteration must be more than 3%, that means the the integration is not reliable). Besides, for a certain sub-process, when some parameters are changed, such as the rapidity range, the p t points, the polarization direction, and so on, the convergence behavior of integration could change at the same time.
In this case, the number of sample points should be reset for the sub-processes to ensure a smooth result in the whole p t region, especially for the loop-diagram sub-processes such as 'g + g →J/ψ + g ', and the sub-process 'g + g →J/ψ + g + g ', since the value of these two processes are usually large and opposite, which could lead to large number canceling. Fig.4 and Fig.5 is a complete result for one p t point. If more p t points were listed in input.dat file, the results will be written into one convergence.dat in the order of p t points, as well as the file fresult.dat.

VIII. SUMMARY
We have presented the Fortran package FDCHQHP, which is used to do the calculation on the p t distribution of yield and polarization for heavy quarkonium prompt production in hadron colliders. The package contains all the sub-processes of CS and CO states up to NLO level, and it focus on the calculation of the polarization parameters.
This package has been fully developed in our past works with rather reliable performance and accuracy. We also provide a small package WBIN, which can help users to run FD-CHQHP in supercomputers easily. Some other improvement are in progress so that we can make it more automatically.