Dealing with redundancy and inconsistency in constructive geometric constraint solving

Available online at www.sciencedirect.com

Advances in Engineering Software 39 (2008) 770–786 www.elsevier.com/locate/advengsoft

Dealing with redundancy and inconsistency in constructive geometric constraint solving David Podgorelec *, Borut Zˇalik, Vid Domiter Faculty of Electrical Engineering and Computer Science, University of Maribor, Smetanova 17, SI-2000 Maribor, Slovenia Received 26 April 2006; received in revised form 15 January 2007; accepted 12 October 2007 Available online 3 December 2007

Abstract General constructive geometric constraint solvers are pre-processed by a degree-of-freedom analysis, which enables efficient graph decomposition and recombination. However, all these methods are based on the assumption that structural rigidity automatically assures solvability. In this paper, we show that this assumption fails in numerous, even the most basic, configurations. We introduce several simple but efficient rules aimed to additionally analyse solvability in such cases. Another novelty addresses conditional constraints between three or more geometric parts, rules for their simplification and a redundancy check. All these functionalities are built into our original 2D geometric constraint solver, based on concepts of rigid clusters and constrained-angle (CA) sets.  2007 Elsevier Ltd. All rights reserved. Keywords: Geometric constraint solving; Computer-aided design; Constructive approach; Redundant constraints

1. Introduction Geometric constraint solving (GCS) represents the heart of modern computer-aided design (CAD). A geometric constraint problem (GCP) consisting of a finite set of geometric elements and a finite set of geometric constraints (relationships between the elements), captures a designer’s intent more naturally than traditional geometric modelling techniques. The solution of the GCP is a class of geometric element instantiations in a given coordinate system, such that all the constraints are satisfied [1]. Of course, it is expected to obtain this solution automatically. Geometric constraint solvers can be classified into two main groups: direct and constructive solvers. A direct solver transforms constraints into equations and, simultaneously, solves them by employing a general numerical [2,3] or symbolic (algebraic) technique [4,5]. Several algorithms were proposed for decomposing sparse polynomial systems of *

Corresponding author. Tel.: +386 2 220 74 74; fax: +386 2 220 72 72. E-mail addresses: [email protected] (D. Podgorelec), zalik@ ˇ alik), [email protected] (V. Domiter). uni-mb.si (B. Z URL: http://gemma.uni-mb.si (D. Podgorelec). 0965-9978/$ - see front matter  2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.advengsoft.2007.10.003

equations [6], but the time complexity remains at least exponential in the size of the largest, independently solved, subproblem. The use of direct solvers should be, therefore, restricted to subproblems that are as small as possible. Isolation of such subproblems, whose solutions can be recombined by solving other small subproblems, represents the main task of the constructive solvers. These solvers typically represent a GCP by a constraint graph G = (N, E) where the set of nodes N consists of geometric elements and the edges from E represent constraints. Local propagation [7] handles each constraint independently. The constraints that should be simultaneously solved form so-called cycles, and these make the problem insoluble. More advanced shape-recognition-based (SR-based) solvers [8–11] use a set of construction rules to solve some predefined subgraph patterns, but they also have limited scope. We can extend it by increasing the repertoire of patterns, but doing so results in greater combinatorial complexity and makes an efficient implementation difficult [1]. Present-day research in GCS tends towards general geometric constraint solvers. These methods are based on a degree-of-freedom (DOF) analysis. In [12], Hoffmann et al. analysed the most prominent representatives of the

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

SR-based and so-called maximum-matching-based (MMbased) solvers. The latter have several advantages, including generality, but they still suffer from some drawbacks. In [13], the same team of authors presented a networkflow-based algorithm for finding minimal dense subgraphs. They employed it in the so-called modified frontier algorithm (MFA) [14]. Gao and Zhang [15] provided another powerful solution with their c-tree decomposition algorithm. They combine the idea of s-tree by Joan-Arinyo et al. [10] and Latham-Middletich’s connectivity analysis [16]. The MFA and the c-tree method can handle 2D and 3D problems. However, the ‘general’ solvers are also not free of drawbacks. They are all based on the assumption that structurally well-constrained, over-constrained and under-constrained problems are, in general, geometrically well-constrained, over-constrained and under-constrained, respectively. Gao and Zhang [15] rightly admit that the DOF analysis actually concerns structural rigidity only. Researchers in GCS are certainly aware of drawbacks and limitations of graph-based methods but, in spite of this, the problem has been neglected and rarely addressed in literature since recently. GCS is closely related to combinatorial rigidity studies, which provide important definitions of rigidity and its applicable approximations, and to automatic geometric theorem proving. Recent works of Sitharam and Zhou [17], Foufou and Michelucci [18,19] successfully link discoveries from the three related fields and significantly extend generality of GCS. However, the automatic geometric theorem proving requires more general techniques that are much more complex than those required by GCS are. Moreover, these techniques typically rely on symbolic or numerical computations that, let us recall, are exactly what we try to avoid with the constructive approach. This means that pure graph-based methods still deserve efforts to improve them to better detect dependences between geometric constraints. Jermann et al. [20] achieved further generalisation of the MFA [13,14] by using an extended definition of structural rigidity, and by distributing flow in a geometrically correct way. Pabon’s idea [21] to classify DOFs as dimensional, translational and rotational also has potential and deserves more attention. In this paper, we introduce an efficient algorithm detecting redundancy, consistency and solvability in 2D GCPs described by a rich, natural, although limited, repertoire of constraints. We are not preoccupied with generality and, therefore, we found it acceptable to use a SR-based solver as a framework for adding new functionalities. Definitions of structural rigidity, geometrical solvability, redundancy and inconsistency are given in Section 2. The structure and the originality of our geometric constraint solver are explained in Section 3. The graph decomposition phase described in Section 4 is a variant of works by Sunde [22], and Verroust et al. [23]. Rigid initial clusters are extracted in this phase, and organised into generally larger constrained-angle sets (CA sets). The recombination phase is considered in Section 5. The set of cluster merging rules

771

is extended in comparison to the repertoires from [8,9,23]. Some special rules handle so-called conditional constraints between three or more geometric elements. Several existing constructive solvers are also capable of handling constraint hypergraphs with each hyperedge involving any number of nodes [12,24]. However, we show here that this topic can also contradict generality. In Section 6, we give several practical examples. Finally, we offer a brief summary in Section 7. 2. Structural rigidity and geometrical solvability Intuitively, rigidity means that presented geometric configuration’s shape and size cannot be changed without modifying its description. Automatic determination of rigidity represents one of the key problems in modern GCS, automatic geometric theorem proving and combinatorial rigidity theory. Since solutions of a GCP do not depend only on a constraint graph structure but also on numerical values of constraints’ parameters, practical applications mostly rely on DOF-based approximation of rigidity called structural rigidity. Definition 1 introduces DOF-based classification of GCPs. Here we assume that each geometric element has two DOFs (like a point in Euclidean plane), and each constraint eliminates a single DOF i.e. it is represented by a single scalar equation (like a distance between two points). The GCP is assumed to retain three DOFs, two translations and a rotation. This makes sense since designers are usually interested in an object’s shape and size, unlike its position and orientation in some fictional coordinate system. Definition 1. The constraint graph G = (N, E) and, simultaneously, the GCP are: • structurally over-constrained if there is an induced subgraph G 0 = (N 0 , E 0 ) with |E 0 | > 2|N 0 |  3, • structurally under-constrained if G is not structurally over-constrained and |E| < 2|N|  3, or • structurally well-constrained if G is not structurally overconstrained and |E| = 2|N|  3.

Definition 2. A constraint e is structurally redundant if we can derive a structurally over-constrained subgraph G 0 = (N 0 , E 0 ) of G, e 2 E 0 , such that G00 = (N 0 , E 0  {e}) is structurally well-constrained. Subgraph G 0 in Definition 2 need not be an induced subgraph of G. Both polygons in Fig. 1 contain structurally redundant constraints. In Fig. 1a, we can derive G 0 = (N 0 , E 0 ), where N 0 = {p1, p2, p3, p5}, and E 0 = {e1, . . ., e6}. In Fig. 1b we can detect N 0 = {p1, p2, p3, p4, p5} and, for example, E 0 = {e1, e2, e3, e4, e7, e8, e9, e11}. In both cases, any constraint ej from E 0 can be considered structurally redundant since G00 = (N 0 , E 0  {ej}) is structurally wellconstrained for all ej. Note that only G 0 in case (a) is an induced subgraph of the corresponding constraint graph G.

772

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

Fig. 1. Two examples of polygons containing structurally redundant constraints.

Structurally redundant constraints must be excluded from the constraint satisfaction and checked for consistency afterwards. A consistently structurally redundant constraint is one which can be added to the constraint system (or removed from it) without changing the set of solutions. Such constraints are automatically satisfied after the complementary subproblem is solved, and there is no need to alarm the user. On the other hand, the inconsistently structurally redundant constraints contradict all the solutions. They should be reported to the user and interactively removed afterwards. Until now, we have considered only structural rigidity. However, the absence of inconsistently structurally redundant constraints does not guarantee that a GCP is solvable. Namely, a geometric configuration can hide various domain-dependent mathematical theorems establishing relations between dimensions, which may be considered additional implicit constraints. Since a user is typically not aware of presence of these implicit constraints, he or she often tries to constrain an internally established relations by additional constraints. We call these constraints geometrically redundant constraints. They are either consistently geometrically redundant or inconsistently geometrically redundant. Detection of geometrical redundancy is the main task of automatic geometric theorem proving. Both, structurally or geometrically redundant constraints, cause that equations, representing a GCP, are not independent. Constraints, that are not redundant, will consequently be called independent constraints. Informally, an over-constrained problem has no solutions, a well-con-

strained problem has a finite number of solutions, but an under-constrained problem has infinitively many solutions [9]. Under the limitations, considered in Definition 1, we are now able to define these terms more formally. Definition 3. A GCP described by the constraint graph G = (N, E) is: • over-constrained (or geometrically over-constrained) if E contains inconsistently (either structurally or geometrically) redundant constraints, • under-constrained (or geometrically under-constrained) if it is not over-constrained and if E contains less then 2|N|  3 independent constraints, • well-constrained (or geometrically well-constrained) if it is not over-constrained and if E contains exactly 2|N|  3 independent constraints, • solvable (or geometrically solvable) if it is well constrained or under-constrained. Fig. 2 shows three structurally well-constrained, but geometrically over-constrained GCPs. Angles a, b, c and d in case: (a) are not independent since the sum of the interior angles in a n-sided polygon is always (n  2)p. In case (b), a = 2b, as a and b represent the central and inscribed angles above the same chord in the circle, respectively. Finally, redundancy in case (c) originates from an internally established relation a2 + c2 = b2 + d2  2ef cos u. Here e and f represent the diagonals. For u = p/2, used in Fig. 2c, the above implicit constraint becomes

Fig. 2. Three examples (a, b, c) of structurally well-constrained, but geometrically over-constrained quadrilaterals.

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

a2 + c2 = b2 + d2. The dimensions a, b, c, d, and u = p/2 are obviously not independent. Dealing with geometrical solvability is evidently much harder and an apparently more common problem than presented so far. One may simply discover additional problematic examples by studying, for example, cyclic or bicentric quadrilaterals. 3. The method The main goal of this research is to extend the drawing scope of our 2D SR-based geometric constraint solver proposed in [25]. That method automatically derives redundant constraints that may be used in local propagation together with or instead of the original ones inserted by a user. Solving triangles, determining the sums of adjacent angles, and elimination of the so-called simple conditional constraints are used to determine redundant distances and angles. The former are stored in the matrix of distances MD, and the latter in the matrix of angles MA. The redundant constraints often ‘break’ cycles, but the drawing scope, although extended in comparison to local propagation, remains limited. Fig. 3 shows two well-constrained quadrilaterals that cannot be solved by the approach from [25]. The case (a) is constructible with ruler and compass, but the solver only manages to determine angle \(l2, l4), which does not prove useful. In the case (b), not even a single redundant constraint is derived. We consider only GCPs in Euclidean plane. The geometric elements are points and lines, which form the so-called auxiliary geometry. We recall that more complex geometric elements of the visible geometry as, for example, circular arcs and cubic Be´zier curves, can be handled by mapping them onto the control points and lines [7,25]. Each point pi is described with its Cartesian coordinates (xi, yi), and each line lj with its slope angle aj 2 [0, p) and a point lying on it. We measure all angles in counter-clockwise direction. The slope angle 0 corresponds to a horizontal line. Algorithm 1 roughly outlines the structure of the presented method. Steps 1 to 4, 10, and 11 have been described in our previous work [25]. Let us briefly describe an acceleration of Step 4, representing the only significant modification in these inherited parts. In the set of n points, we can define O(n2) line segments, and O(n4) angles between

773

the line segments. Obviously, O(n4) time is needed for a single traverse of all the angles. In a brute-force approach, the number of newly determined angles in a single traverse usually significantly decreases after a few iterations, resulting in a catastrophic worst-case time complexity O(n8). On the other hand, the accelerated approach iteratively traverses the angles only a few times, but then switches to another technique propagating the new angles in the graph G 0 = (N 0 , E 0 ). Here the nodes of N 0 are line segments and the edges of V 0 are the angles. In this way, we retain the expected time complexity O(n4), which does not exceed the worst-case time complexity of Algorithm 1. Algorithm 1. Geometric constraint solver 1. Mapping the visible geometry onto the auxiliary geometry. 2. Transformation of composed constraints into simpler ones. 3. Splitting the constraint graph into independently solvable connected subgraphs (subproblems). for (each independently solvable connected subgraph) do begin 4. Adding redundant distances and angles. 5. Transforming the subproblem to the form addressing points’ coordinates only. 6. Creation of initial clusters. 7. Cluster registration. 8. Organising registered clusters into initial CA sets. 9. Recombination phase. 10. Absolute positioning. end for 11. Re-establishment of the visible geometry from the positioned auxiliary geometry. Step 5 is described in Section 3.1. It transforms a rich repertoire of constraints into the equations on points’ coordinates only. Thus, the nodes of the constraint graph used in the subsequent steps are actually the points only. Steps 6, 7, and 8, representing the graph decomposition phase, are explained in Section 4. This phase is based on the works by Sunde [22], and Verroust et al. [23], but pre-processing

Fig. 3. Calculation of redundant dimensions does not provide for the solvability of quadrilaterals (a, b).

774

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

by adding redundant constraints entails some original features. In addition, the cluster registration which importantly reduces the number of clusters for further processing, represents a completely original concept. The recombination phase considered in Section 5 brings several original contributions. We introduce the hierarchical structure, which enables calculation of only those redundant dimensions that are really needed. The set of cluster merging rules is extended in comparison to the repertoires used in [8,9,23]. Special treatment of angular conditional constraints often offer new possibilities for clusters or CA sets merging. Finally, we propose an original heuristic approach for determining the structural rigidity. 3.1. Transformation of constraints into equations

8. Sum4L(l1, l2, l3, l4, a) means \(l1, l2) + \(l3, l4) = a. The equation is:u(p1, p2, p3, p4) u(p5, p6, p7, p8)  (1  u 2 (p 1 , p 2 , p 3 , p 4 )) 1/2 (1  u 2 (p 5 , p 6 , p 7 , p 8 )) 1/2  cos a = 0. 9. Sum6L(l1, l2, l3, l4, l5, l6) means \(l1, l2) + \(l3, l4) = \(l5, l6). The equation is: u(p1, p2, p3, p4) u(p5, p6, p7, p8)  (1  u2(p1, p2, p3, p4))1/2(1  u2(p5, p6, p7, p8))1/2  u(p9, p10, p11, p12) = 0. 10. RatioL(l1, l2, l3, l4, r1, r2) means \(l1, l2)/\(l3, l4) = r1/ r2. The equation is:   br =2c 2 X r2 k k ð1Þ ur2 2k ðp1 ; p2 ; p3 ; p4 Þð1  u2 ðp1 ; p2 ; p3 ; p4 ÞÞ 2k k¼0   br 1 =2c X r1 k  ð1Þ ur1 2k ðp5 ; p6 ; p7 ; p8 Þ 2k k¼0 k

A user can use 26 predefined constraint types, but constraints of 15 types are decomposed into conjunctions of intuitively simpler constraints in Step 2 of Algorithm 1. Only 11 constraint types are therefore used in a transformed GCP. Since we wish the GCP to be solved as simply as possible, Step 5 of Algorithm 1 transforms constraints addressing lines into constraints among the corresponding control points. In continuation, we give the equations on points’ coordinates, describing the transformed 11 constraint types. Each line lj in constraints 3, 8, 9, or 10 passes through a pair of points (p2j 1, p2j). Presentation of the distance between a pair of points is based on the formula for the Euclidean distance in plane: d(p1, p2) = ((x1  x2)2 + (y1  y2)2)1/2. Transformation of constraints addressing angles between lines is based on the scalar product. If l1 = (p1, p2), and l2 = (p3, p4), then function u(p1, p2, p3, p4) = ((x2  x1)(x4  x3) + (y2  y1)(y4  y3))/(|p1p2||p3p4|) equals cos \(l1, l2). The function is applicable when |p1p2| and |p3p4| are both known. If we do not know a line segment with known length along any of the lines involved in the relationship, then the method inserts a new point of the so-called invisible geometry, and places it a unit away from some known point on the line. Such a point remains invisible to the user and, therefore, he/she cannot use it in other constraints. More details on invisible geometry can be read in [25]. 1. Distance(p1, p2, d), d > 0 (if d = 0 then Coincidence is used), is represented by: d(p1, p2)  d = 0. 2. Coincidence(p1, p2) is the only constraint requiring two equations: x1  x2 = 0, y1  y2 = 0. 3. Angle(l1, l2, a) introduces the relationship \(l1, l2) = a. The equation is: u(p1, p2, p3, p4)  cos a = 0. 4. Difference(p1, p2, p3, p4, d) means |p1p2|  |p3p4| = d. The equation is: d(p1, p2)  d(p3, p4)  d = 0. 5. Sum4(p1, p2, p3, p4, s) means |p1p2| + |p3p4| = s. The equation is: d(p1, p2) + d(p3, p4)  s = 0. 6. Sum6(p1, p2, p3, p4, p5, p6) means |p1p2| + |p3p4| = |p5p6| i.e. d(p1, p2) + d(p3, p4)  d(p5, p6) = 0. 7. Ratio(p1, p2, p3, p4, r) means |p1p2|/|p3p4| = r. The equation is: d(p1, p2)  rd(p3, p4) = 0.

 ð1  u2 ðp5 ; p6 ; p7 ; p8 ÞÞ ¼ 0: 11. Collinear(p1, p2, p3) actually represents Sum6(p1, p2, p2, p3, p1, p3). However, the overall treatment of collinearity involves some particularities that are uncharacteristic for other appearances of Sum6. The constraints 4–11 are conditional constraints. Each of them may be replaced by a simple dimensional constraint Distance or Angle if all other dimensions in the relationship are known. If we introduce a uniform numbering of variables, for example by mapping yk ! xm+k where m is the number of points in the observed equation, then all the equations PQ nðj;iÞ 4 are in form ðxÞ ¼ 0, where x represents the j¼1 fj;i i vector of variables (x1, . . ., x2m), the exponents n(j, i) are rational numbers and all the functions fj,i(x) are algebraic. Furthermore, this form is also applicable for all partial derivatives of the first and the second orders, usually required by employed numerical methods [3]. A unified representation of functions and derivatives in a rather simple form results in at least two advantages: (1) all partial derivatives can be simply computed in a symbolic way by employing a special data structure, and (2) development of an efficient user interface will enable users to interactively define new constraint types. 4. Decomposition phase The graph decomposition algorithm presented by Sunde [22], and later also employed by Verroust et al. [23] operates with two types of sets of geometric elements. A constrained angle set (in short, CA set) is a set of pairs of points whose corresponding oriented segments are mutually constrained in angle. A constrained distance set (CD set) is a set of points with mutually constrained distances. Elementary CA and CD sets are created when constraints are added to the model. When an angle constraint between two directed segments is inserted, a CA set is created, and the orientations are noted. When a distance constraint between two points is given, a CD set and a CA set are cre-

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

ated. These elementary CA and CD sets are later merged using the construction rules [23]. Our graph decomposition algorithm represents a variant of the described method. We define the incorporated datasets in a modified, more general manner. We replace CD sets by clusters but we keep labelling CD1, CD2, . . ., CDm for the clusters. We also retain CA sets. Definition 4. A cluster is a maximal structurally rigid geometric part. The required maximality in Definition 4 means that none of clusters may represent a subset of any other cluster. The requirement for structural rigidity implies two useful corollaries. Corollary 1. Any distance or angle inside a cluster is known. It is either explicitly defined or computable in O(1) time from other known dimensions in the cluster. Corollary 2. A cluster is easily constructible by local propagation in linear time. Definition 4 of clusters is much more flexible than Sunde’s definition of CD sets. We may say that a cluster contains points, line segments, pairs of line segments, triangles etc. With a CD set, any generalisation of this kind is formally unacceptable because the CD sets are defined as point sets [22,23]. There is at least one more difference between clusters and CD sets. We require that each point of the auxiliary and the invisible geometry should be a member of at least one cluster. Since a GCP may contain points that are not addressed by any distance constraints, a cluster containing a single point is also acceptable. We recall that an elementary CD set contains two points. Definition 5. A cluster is called 1-cluster if it contains a single point, 2-cluster if it contains two points, 3-cluster if it contains three points, or k-cluster if it contains three or more points. 1-clusters usually require special treatment since they have only 2 DOFs instead of 3. For brevity, we will mostly deal only with 2-clusters and k-clusters. Note that each 3cluster is the k-cluster as well. Different choices of initial clusters are possible. We assume that any distance or angle inside an initial cluster is explicitly defined before the decomposition phase starts. Since the decomposition is pre-processed by adding redundant dimensions [25], each initial cluster is either a CD set or 1-cluster. Theorem 1. Two initial clusters can share at most one common point. With other words, a line segment with known length cannot belong to two different initial clusters. Theorem 2. If two points pi and ps belong to some initial cluster (CDr) and if there is the third point pk, such that |pipk| and |pspk| are known, then pk should also be inserted into CDr.

775

Theorems 1 and 2 describe two trivial consequences of the mentioned pre-processing, that prove extremely useful in Algorithm 2 responsible for creation of initial clusters. The algorithm assumes that the distances are sorted with respect to the points’ indices. The use of the matrix of distances MD established in the pre-processing phase [25] guarantees the required presortedness. Functions first(CDr) and second(CDr) return points with the lowest and the second lowest index in CDr, respectively. Algorithm 2. Cluster Creation. Initialise empty lists of clusters. Label all the known distances in the matrix of distances MD as NOT_USED. for (each point pi) do begin CD :¼ {pi, pj}, where |pipj|, j > i, is the first distance (with the lowest j), labelled as NOT_USED. if (CD established) then begin Add CD to the list of clusters, and label |pipj| as USED. for (all |pipk|, k > j, labelled as NOT_USED) do begin for (each cluster CDr, where first (CDr) = pi) do begin ps :¼ second(CDr). if (|pspk| is known) then begin (* use of Theorem 2 *) CDr :¼ CDr [ {pk}. Label the distances between pk and all other points from CDr as USED. exit for. end if end for. if (pk has not been added to any cluster) then Add {pi, pk} to the list of clusters, and label |pipk| as USED. end for end if. if (pi does not belong to any cluster) then Add {pi} to the list of clusters. end for.

4.1. Cluster registration Definition 6. A constraint is called an inter-cluster constraint if it addresses points from different clusters. Constraints that are not inter-cluster will be called the internal cluster constraints. Algorithm 2 treats all the distances as internal cluster constraints. Consequently, it usually produces a high number of initial 2-clusters. For example, the well-known smallest tri-connected graph in Fig. 4(a) is decomposed into three 2-clusters and two k-clusters. The example in

776

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

Fig. 4. Cluster registration reduces the number of clusters from 5 to 2 (a, b), and from 32 to 8 (c, d).

Fig. 4c is even more expressive. Each of 32 constraints Distance results in a separate 2-cluster. We will show, in continuation, that the recombination phase analytically solves only patterns with up to three clusters and, therefore, the obtained decomposition is not optimal for further processing. We introduce the cluster registration algorithm to prevent such situations. It registers all 1-clusters and k-clusters, and those 2-clusters, assuring that each point remains a member of at least one cluster. Only these registered clusters enter the recombination phase, and all other 2-clusters are replaced by inter-cluster constraints Distance. Fig. 4b and d show only the registered clusters, corresponding to the problems from Fig. 4a and c, respectively. Algorithm 3. Cluster registration. 1. Register all 1-clusters and k-clusters. 2. for (each 2-cluster CDi with both points found in registered k-clusters) do 3. Replace CDi by the corresponding inter-cluster constraint Distance. 4. Let P represent the set of points that do not belong to any registered cluster. 5. while (P is not empty) do begin 6. Select a point pi 2 P, and remove it from P. 7. Register a 2-cluster containing pi and some other point pj.

8. if (pj 2 P) then Remove pj from P. 9. end while. 10. for (each 2-cluster CDi) do 11. if ((CDi is not registered) or (CDi shares its two points with a pair of registered clusters)) then 12. Replace CDi by the corresponding intercluster constraint Distance. Note that Fig. 4d does not represent the only possible solution. Eight nearly horizontal or eight non-intersecting diagonal 2-clusters are also acceptable, and there are many other solutions with eight or more registered clusters. The result depends on initial clusters indexing, and on rules for selecting pi and pj in clauses 6 and 7 of Algorithm 3. However, no matter which criteria we use for these selections, Algorithm 3 assures that each registered 2cluster contains at least one point which is not present in any other registered cluster. This enables us to state the extremely important Theorem 3. Theorem 3. The number of registered clusters cannot exceed the number of points n. Proof. Let us assume the opposite: a GCP with n points is decomposed into m registered clusters, where 0 < n < m. Let r 6 m represent the total number of 1-clusters and 2-clusters. If we remove these r registered clusters, we

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

simultaneously remove at least r points. We obtain a subproblem with, at most, p = n  r points and with exactly q = m  r k-clusters. Our hypothesis 0 < n < m becomes 0 6 p < q. Then we reduce each k-cluster into a 3-cluster. This operation may also reduce the point set, but the number of 3-clusters (former k-clusters) remains q. Let N 0 and E 0 represent the point set and the set of the internal cluster constraints Distance in all 3-clusters, respectively. There are three such constraints in each 3-cluster and, consequently, |E 0 | = 3q. GCP 0 = (N 0 , E 0 ) is certainly not structurally over-constrained because none of the triangles (3-clusters) constrained by three side lengths is not structurally over-constrained. Therefore, |E 0 | 6 2|N 0 |  3. We now use 0 6 |N 0 | 6 p < q, derived from the initial hypothesis, and |E 0 | = 3q. We obtain 3q < 2q  3. The solution q < 3 proves that the initial hypothesis was incorrect since the number of k-clusters is not allowed to be negative. h

777

Create a new CA set CAi, insert CDr into CAi, and label CDr as USED. for (each cluster CDj with known angle between it and CDr) do Insert CDj into CAi, and label CDj as USED. forever.

Theorem 5. Two clusters belonging to the same initial CA set cannot share a point. Theorem 6. Each vector, each initial 2-cluster or k-cluster belongs to exactly one initial CA set. Theorem 7. Initial CA sets are maximal i.e. none of them represents a subset of any other initial CA set. 5. Recombination phase

4.2. Organising clusters into CA sets Here we define CA sets in more general way than Sunde [22] and Verroust et al. [23] did. Then we give a straightforward algorithm for creation of initial CA sets. Our definition is based on terms of vectors and angles between clusters that need be introduced first. Definition 7. A pair of points with known mutual distance will be called a vector. Theorem 4. When the angle between a pair of vectors vi 2 CDq and v0i 2 CDr is constrained, the angles between all the pairs of vectors vj 2 CDq and v0j 2 CDr are constrained. The angle \ðvi ; v0i Þ may be simply referred to as the angle between the two clusters. Proof. \ðvj ; v0j Þ ¼ \ðvj ; vi Þ þ \ðvi ; v0i Þ þ \ðv0i ; v0j Þ. The first and the last summand on the right side are known since they represent internal cluster constraints. The sum is therefore known if the second summand is also known. h Definition 8. A constrained angle set (CA set) is a 2-cluster, a k-cluster or a structurally under-constrained geometric part consisting of 2-clusters and/or k-clusters with known angles between them. Definition 8 gives straightforward instructions how to design Algorithm 4 for CA sets creation. Theorems 5–7 stress some properties of initial CA sets. The proofs are trivial and left to the reader.

Algorithm 4. Creation of CA Sets. Label all 2-clusters and k-clusters as NOT_USED. loop Select a cluster CDr, labelled as NOT_USED. if (nothing selected) then exit loop.

The recombination phase merges clusters into larger clusters, and CA sets into larger CA sets. A model is completely constrained when all the points belong to the same cluster. The algorithm searches for pairs or triples of clusters that can be merged by employing some pre-defined construction rules. It traverses all O(m2) pairs of clusters and all the triples containing some previously identified adjacent pairs. This process can be repeated at most m times, resulting in the worst-case time complexity O(m4). Theorem 3 assures that the number of clusters m belongs to O(n) but, typically, m is much smaller than the number of points n. Simultaneously with the cluster merging, the process of merging CA sets is also running. This process usually results in newly determined angles between pairs of clusters from the merged CA sets, and these angles may then encourage further cluster merging. If the described process does not result in a single cluster then we must cope with larger k-tuples of clusters. The brute force approach traverses all possible combinations. This strategy is reliable but it requires exponential time and is, therefore, applicable only when the number of clusters is considerably low. For this reason, we introduce a more elegant, much quicker but less reliable heuristic approach. Algorithm 5 outlines the structure of the recombination phase. In a practical implementation, there is much higher coherence between the particular steps, but these details are beyond the scope of this paper. Sections 5.1– 5.5 describe particular steps of the algorithm, the supporting data structure and its management. In Definitions 9 and 10, we define some terms, basic for understanding these Sections. Definition 9. A point is called a frontier point if it appears in any inter-cluster constraint or is shared by distinct clusters. The active frontier points are those used in the observed merging operation.

778

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

Definition 10. Clusters CD1, CD2, . . ., CDm (or CA sets CA1, . . ., CAk) are linked by constraint ei when all the frontier points addressed by ei belong to CD1 [ CD2 [ . . . [ CDm (or CA1 [ . . . [ CAk). Algorithm 5. Recombination Phase. loop loop 1. Insert redundant dimensions. 2. Perform construction rules to merge clusters and CA sets. if (a single cluster is obtained or no rule has been employed in the last iteration) then exit loop. forever. if (a single cluster is obtained) then exit loop. 3. Run the heuristic algorithm to search for a minimal mergeable subproblem. 4. Numerically construct the minimal mergeable subproblem. if (a single cluster is obtained) then exit loop. forever.

5.1. Calculation of redundant dimensions Let us observe the inner loop of Algorithm 5. Cluster merging and CA set merging operations usually enable us to determine additional redundant distances and angles that have not been identified in the pre-processing phase (step 4 of Algorithm 1). These redundant dimensions then often launch further merging operations. Let us recall from [25] that solving triangles, determining the sums of adjacent angles, and elimination of the simple conditional constraints are used to determine redundant distances and angles. Here we introduce some additional simple rules for the transformation, or even elimination, of conditional constraints addressing the frontier points from exactly two CA sets. It is always worth replacing conditional constraints by simple dimensional constraints whenever possi-

ble. A distance or an angle always link at most two clusters. On the other hand, the current constraint set contains conditional constraints that may link up to six clusters. Theorem 8 generalises the use of Theorem 4 from clusters to CA sets. It provides a partial result for the much more important theorems 9 to 11. Vectors v1, v2 and v3 are from a CA set CA1, and vectors v01 , v02 and v03 belong to another CA set CA2. Theorem 8. The difference between the undetermined angles \ðv1 ; v01 Þ and \ðv2 ; v02 Þ is known. It is equal to the difference between the pair of known angles: \ðv1 ; v01 Þ  \ðv2 ; v02 Þ ¼ \ðv1 ; v2 Þ  \ðv01 ; v02 Þ. Proof. We may simply break \(v1,v 0 1) into the sum: \ðv1 ; v01 Þ ¼ \ðv1 ; v2 Þ þ \ðv2 ; v02 Þ þ \ðv02 ; v01 Þ. The desired relationship is then obtained by subtracting \ðv2 ; v02 Þ on both sides of the equation, and by taking into the account \ðv02 ; v01 Þ ¼ \ðv01 ; v02 Þ. h Theorem 9. Constraint Sum4L(l1, l2, l3, l4, a) linking two CA sets is either (a) redundant or (b) may be replaced by some newly computed angle. Proof. (a) Let v1 2 l1, v01 2 l2 , v02 2 l3 , and v2 2 l4. Sum4L(l1, l2, l3, l4, a) then establishes the relation \ðv1 ; v01 Þ þ \ðv02 ; v2 Þ ¼ a. This is equivalent to \ðv1 ; v01 Þ \ðv2 ; v02 Þ ¼ a. Theorem 8 assures the left part of this equation is already known and, therefore, the constraint is apparently redundant.(b) Let v1 2 l1, v01 2 l2 , v2 2 l3, and v02 2 l4 . Sum4L(l1, l2, l3, l4, a) then establishes the relation \ðv1 ; v01 Þ þ \ðv2 ; v02 Þ ¼ a. We may use Theorem 8 to get rid of \ðv2 ; v02 Þ. The resulting newly computed angle is \ðv1 ; v01 Þ ¼ ð\ðv1 , v2 Þ  \ðv01 ; v02 Þ þ aÞ=2. h Theorem 10. Constraint Sum6L(l1, l2, l3, l4, l5, l6) linking two CA sets may be replaced by the angle. Proof. We must consider four distinct cases. For brevity, we only classify them and give the resulting angles. Arrows in Fig. 5 point from the first vector forming the observed angle towards the second one.

Fig. 5. Theorem 10 requires separate consideration of four different cases (a–d).

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

(a) Let v1 2 l1, v01 2 l2 , v2 2 l3, v02 2 l4 , v3 2 l5, and v03 2 l6 , as shown in Fig. 5a. The equation is \ðv1 ; v01 Þ þ \ðv2 ; v02 Þ ¼ \ðv3 ; v03 Þ, and the resulting new angle is \ðv1 ; v01 Þ ¼ \ðv3 ; v2 Þ  \ðv03 ; v02 Þ. (b) Let v1 2 l1, v01 2 l2 ; v2 2 l3 ; v02 2 l4 ; v03 2 l5 , and v3 2 l6, as shown in Fig. 5b. The constraint equation is \ðv1 ; v01 Þ þ \ðv2 ; v02 Þ ¼ \ðv03 ; v3 Þ, and the resulting new angle is \ðv1 ; v01 Þ ¼ ð\ðv2 ; v1 Þ  \ðv02 ; v01 Þ 0 0 þ\ðv3 ; v1 Þ  \ðv3 ; v1 ÞÞ=3. (c) Let v1 2 l1, v01 2 l2 , v02 2 l3 ; v2 2 l4 ; v3 2 l5 , and v03 2 l6 , as shown in Fig. 5c. The equation is \ðv1 ; v01 Þ þ\ðv02 ; v2 Þ ¼ \ðv3 ; v03 Þ, and the resulting new angle is \ðv3 ; v03 Þ ¼ \ðv1 ; v2 Þ  \ðv01 ; v02 Þ. (d) Let v1 2 l1, v01 2 l2 , v02 2 l3 , v2 2 l4, v03 2 l5 , and v3 2 l6, as shown in Fig. 5d. The equation is \ðv1 ; v01 Þ þ\ðv02 ; v2 Þ ¼ \ðv03 ; v3 Þ, and the resulting new angle is \ðv03 ; v3 Þ ¼ \ðv1 ; v2 Þ  \ðv01 ; v02 Þ. h Theorem 11. Constraint RatioL(l1, l2, l3, l4, r1, r2) linking two CA sets may be replaced by the angle. Proof. (a) Let v1 2 l1 ; v01 2 l2 , v2 2 l3, and v02 2 l4 . The equation is \ðv1 ; v01 Þ=\ðv2 ; v02 Þ ¼ r1 =r2 , and the resulting newly determined angle is \ðv1 ; v01 Þ ¼ r1 ð\ðv1 ; v2 Þ \ðv01 ; v02 ÞÞ=ðr1  r2 Þ. (b) Let v1 2 l1, v01 2 l2 , v02 2 l3 , and v2 2 l4. The equation is \ðv1 ; v01 Þ=\ðv02 ; v2 Þ ¼ r1 =r2 , and the resulting newly determined angle is \ðv1 ; v01 Þ ¼ r1 ð\ðv1 ; v2 Þ 0 0 \ðv1 ; v2 ÞÞ=ðr2  r1 Þ. h Theorem 12. If two clusters CDi and CDj contain a pair of collinear vectors vh 2 CDi and vk 2 CDj then one of the intercluster constraints Collinear may be replaced by the angle between CDi and CDj. Theorems 9–12 offer clear instructions for the determination of new angles. Each of these newly determined angles links two distinct CA sets and launches their merging. The proof of Theorem 12 is trivial. The angle between the collinear vectors is always 0 or p (the appropriate one is chosen by considering the sketch).

5.2. Construction rules Here we describe the repertoire of the employed construction rules for merging clusters and CA sets (step 2 of Algorithm 5). There are three rules aimed to merge pairs of clusters, four rules to merge triples of clusters, and two rules for merging pairs of CA sets. We also use a special rule (Rule 8) that merges a pair of non-rigid geometric parts under certain conditions. Three different techniques of constraint satisfaction are used in our construction rules: local propagation, analytical approach to solving the system of equations, and numerical solving of the system of equations. We favour the use of local propagation because of its simplicity and efficiency.

779

Rule 1: We merge two clusters when they share two (or more) points. See Fig. 6a. This rule cannot be used on initial clusters, but it often becomes useful after employing some other rule. Local propagation is powerful enough to handle these types of problems. Rule 2: We merge two clusters when they share a point and they are linked by a constraint (Fig. 6b). We have provided analytical solutions for this rule by transforming the corresponding systems of equations into a triangular form. Rule 3: We merge two clusters linked by three independent inter-cluster constraints (Fig. 6c). We have provided analytical solutions for those cases when e1, e2, and e3 represent distances, or one of them is an angle and the remaining two are distances. In all other cases, the problem is solved numerically. Rule 4: Let CD1\ CD2 = {pi}, CD2 \ CD3 = {pj}, CD1 \ CD3 = {pk}, i 5 j, j 5 k, i 5 k. We may then merge CD1, CD2, and CD3 (Fig. 6e). Rule 4 represents the basic construction in [9]. Local propagation can handle these types of problems. Rule 5: Let CD1 \ CD2 = {pi}, CD1 \ CD3 = {pj}, i 5 j, and CD2 \ CD3 = ;. If CD2 and CD3 are linked by two independent inter-cluster constraints then we merge CD1, CD2 and CD3 (Fig. 6f). The use of Rule 5 is illustrated in Fig. 6f. Firstly, we merge CD2 and CD3 by employing Rule 3. Here the internal cluster constraint e3 = |pipj| from CD1 substitutes the missing inter-cluster constraint linking CD2 and CD3. After this, Rule 1 is employed to merge (CD2 [ CD3) and CD1. Rule 6: Let CD1 \ CD2 = {pi}, CD1 \ CD3 = ;, and CD2 \ CD3 = ;. If CD1 and CD3 are linked by two independent inter-cluster constraints, and if CD2 and CD3 are linked by two independent inter-cluster constraints then we merge CD1, CD2 and CD3 (Fig. 6g). Rule 7: If six independent inter-cluster constraints addresses CD1, CD2, and CD3, such that each pair of clusters is linked by a pair of constraints, then we merge CD1, CD2 and CD3 (Fig. 6h). In the current implementation, Rule 6 and Rule 7 are handled numerically. Rule 8: If CDi contains all the frontier points adjacent to CDj or shared by CDj, then we merge CDi and CDj (Fig. 6d). Rule 8 represents an exception since it handles structurally under-constrained subproblems. The resulting sub-Figure is not rigid, but it may be later considered and processed as a single cluster because it is not addressed by any inter-cluster constraints. Addition of an inter-cluster constraint therefore requires a consistency test for such ‘‘pseudocluster’’ while, with regular clusters, a constraint deletion is the only critical operation. Fig. 6d shows three possible cases: (1) CD1 and CD4 share point p1, (2) CD2 and CD5 are linked by an inter-cluster constraint e1, and (3) CD3 and CD5 are linked by two inter-cluster constraints e2 and e3. In the first case, only a translation of one of the clusters is required. The second case is handled analytically, regardless of the type of e1. The last case is handled analytically if e2 and e3 represent distances, or

780

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

Fig. 6. Merging (a–d) two or (e–h) three clusters by employing Rules 1–8.

a distance and an angle. Otherwise, the problem is solved numerically. Rule 9: If two CA sets are linked by a constraint Angle then we merge them into a larger CA set. This rule emphasises the importance of Theorems 9–12 and other ways for angles calculation. Rule 10: When clusters from different CA sets are merged, the enclosing CA sets are automatically merged. More than two CA sets can be handled by sequentially adding each of them to the merged CA set. Step 2 of Algorithm 5 tries to use Rule 9 on each pair of CA sets before it employs Rules 1–8 on pairs or triples of clusters. On the other hand, Rule 10 is automatically performed in this second part. It prevents waiting for another iteration of the inner loop of Algorithm 5 and, usually, it slightly accelerates the merging process. Rule 10 is also used during the construction of a minimal mergeable subproblem (step 4 of Algorithm 5) detected by the heuristic method. 5.3. Heuristic search for a minimal mergeable subproblem Before we perform a cluster merging operation, we must be sure that the observed subproblem does not contain smaller mergeable subproblems. Rules 1–8 are designed for patterns that automatically fulfil this requirement but, generally, we have to search for a so-called minimal mergeable subproblem.

Definition 11. A structurally over-constrained subproblem of GCP is called structurally well over-constrained if it can be transformed into a structurally well-constrained one by removing structurally redundant constraints. Definition 12. A subproblem is called minimal mergeable subproblem if it is either structurally well-constrained or structurally well over-constrained, and it does not contain smaller structurally well-constrained or structurally overconstrained subproblems. If we assume that the graph is rid of geometrical redundancy then we could handle this problem by some general constructive solver as, for example, MFA or the c-tree method [14,15]. Here we present a simple heuristic approach that also serves quite well and remarkably fast for this purpose. It realises step 3 of Algorithm 5. The algorithm sequentially eliminates clusters from the graph and checks the structural rigidity of the remaining subgraph. If a structurally well-constrained or a structurally well over-constrained subgraph is identified then we extract it and deal with it independently. Let us suppose that a subgraph contains m active frontier points organised into k clusters. If it is structurally well-constrained or structurally over-constrained then it must contain at least 2m  3 constraints. 2m  3k internal cluster constraints provide structural rigidity of the clusters, and the surplus of 3(k  1) constraints represents the inter-cluster constraints. This value is updated and compared to the number of scalar

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

equations NSE after each cluster elimination. NSE nearly corresponds to the number of inter-cluster constraints, but there are three situations that require special treatment: (1) the constraint Coincidence contributes 2 to NSE, (2) a point shared by k clusters contributes 2(k  1) to NSE, and (3) in a group of k clusters from the same CA set, k(k  1)/ 2 angles between pairs of clusters contribute k  1 to NSE (see Section 5.5 for explanation). The heuristic criterion for selection of a cluster for elimination is based on the clusters’ valences. The valence generally means the number of scalar equations addressing the cluster’s active frontier points, but here we also meet two exceptions: (1) a point shared by k clusters contributes 2(k  1)/k to the valence of each of the k clusters, and (2) in a group of k clusters from the same CA set, angles between these clusters increase valences of each of the k clusters for (k  1)/k when k > 2, or for 1 when k = 2 (see Section 5.5). Let us remark that we may ignore those clusters with the valence 3 or less. Three equations are only sufficient to position a single cluster and cannot contribute to positioning of the remaining clusters. The heuristic criterion consists of six rules listed below. If the first rule finds more candidates for elimination then the second rule is employed etc. The priorities of the rules have been experimentally established during some hundreds of tests. 1. The cluster with the lowest valence in the updated subgraph has to be eliminated. 2. The cluster with the lowest number of neighbours in the updated subgraph has to be eliminated. 3. The cluster updated the most recently has to be eliminated. 4. The cluster with the lowest initial valence has to be eliminated. 5. The cluster with the lowest initial number of neighbours has to be eliminated. 6. The cluster with the lowest index has to be eliminated.

781

The goal of cluster elimination is to isolate a minimal mergeable subproblem by removing structurally underconstrained subproblems. The heuristic criterion is based on the assumption that a structurally under-constrained subgraph is intuitively less dense than the minimal mergeable subgraph. A cluster with low valence and a low number of neighbours is, therefore, a good candidate for belonging to some structurally under-constrained part. The third rule is particularly interesting. It tries to hold the elimination process in a presumably structurally under-constrained part close to the preceding eliminations, instead of jumping chaotically from one end of the graph to another. However, the experimentally obtained heuristic criterion does not guarantee success in general. For this reason, we still separately consider patterns of two and three clusters as described in Section 5.2. Let us consider the three examples presented in Fig. 7a– c. Each of them consists of 10 clusters and 27 inter-cluster constraints. Each triple in the sequence of clusters in Tables 1–3 represents the valence, the number of neighbours, and the cluster’s index. Example. In the graph from Fig. 7a, the process terminates when the under-constrained subproblem with two clusters and NSE = 2 is obtained, as shown in Table 1. NSE was less than 3(k  1) all the time, except at the beginning. The entire problem, therefore, represents the minimal mergeable subproblem. Example. In the graph from Fig. 7(b), the process determines equality NSE = 3(k  1) = 12 for the subproblem with k = 5 clusters, as shown in Table 2. The elimination process proceeds and determines that this subproblem does not contain smaller mergeable subproblems. The further steps are: 1. We numerically merge the structurally well-constrained subproblem consisting of CD6 to CD10.

Fig. 7. A dense (a) and a decomposoble (b) structurally well-constrained problem, a structurally over-constrained problem (c), and a partial result during reduction of the previous problem (d).

782

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

Table 1 Sequence of cluster eliminations in analysing the graph from Fig. 7a Elimination

k

NSE

3(k  1)

Sorted sequence of clusters

– CD2 CD10 CD6 CD1 CD3 CD4, CD5, CD9

10 9 8 7 6 5 2

27 23 20 16 13 9 2

27 24 21 18 15 12 3

(4|3|2), (4|3|4), (4|3|8), (4|3|10), (6|4|1), (6|4|5), (6|4|6), (6|4|7), (7|5|3), (7|5|9) (3|2|10), (4|3|1), (4|3|4), (4|3|8), (6|4|3), (6|4|5), (4|3|6), (6|4|7), (6|4|9) (4|3|6), (4|3|1), (4|3|4), (4|3|8), (6|4|9), (6|4|3), (6|4|5), (6|4|7) (3|2|1), (4|3|9), (4|3|4), (4|3|8), (5|3|7), (6|4|3), (6|4|5) (4|3|3), (4|3|9), (4|3|4), (4|3|8), (5|3|5), (5|3|7) (3|2|4), (3|2|5), (3|2|9), (4|3|8), (5|3|7) (2|1|8), (2|1|7)

Table 2 Sequence of cluster eliminations in analysing the graph from Fig. 7b Elimination

k

NSE

3(k  1)

Sorted sequence of clusters

– CD5 CD2 CD1 CD3, CD4 CD10 CD7, CD9

10 9 8 7 5 4 2

27 23 19 15 12 8 2

27 24 21 18 12 9 3

(4|3|5), (4|3|7), (5|4|2), (5|4|10), (6|4|1), (4|3|2), (4|3|4), (4|3|7), (5|3|1), (5|4|10), (4|2|1), (4|3|10), (4|3|3), (4|3|4), (4|3|7), (2|2|3), (2|2|4), (4|3|10), (4|3|7), (6|4|6), (4|3|10), (4|3|7), (5|3|8), (5|3|9), (6|4|6) (3|2|7), (3|2|9), (5|3|8), (5|3|6) (2|1|6), (2|1|8)

(6|4|3), (6|4|3), (6|4|6), (6|4|8),

(6|4|4), (6|4|6), (6|4|8), (6|4|9) (6|4|6), (6|4|8), (6|4|9) (6|4|8), (6|4|9) (6|4|9)

(6|4|1), (5|4|2), (5|4|2), (5|4|2), (6|4|3)

(6|4|3), (6|5|2), (7|5|4), (7|5|9) (6|4|3), (5|4|4), (7|5|9) (6|4|9), (6|4|3) (6|4|3)

Table 3 Sequence of cluster eliminations in analysing the graph from Fig. 7c Elimination

k

NSE

3(k  1)

Sorted sequence of clusters

– CD5 CD8 CD7 CD6 CD10 CD9 CD2 CD3

10 9 8 7 6 5 4 3 2

27 23 19 16 12 10 9 5 2

27 24 21 18 15 12 9 6 3

(4|3|5), (4|3|8), (4|3|10), (5|3|6), (5|3|7), (4|3|8), (4|3|10), (5|3|1), (5|3|6), (5|3|7), (3|2|7), (4|3|4), (4|3|10), (5|3|1), (5|3|6), (4|2|6), (4|3|9), (4|3|4), (4|3|10), (5|3|1), (2|2|10), (2|2|9), (4|3|4), (5|3|1), (5|4|2), (1|1|9), (4|3|2), (4|3|4), (5|3|1), (6|4|3) (4|3|2), (4|3|4), (5|3|3), (5|3|1) (3|2|3), (3|2|4), (4|2|1) (2|1|1), (2|1|4)

2. The merged cluster CD6to10 may be ignored during the analysis of the simplified graph since its valence is 3. We, therefore, deal with the subproblem consisting of CD1 to CD5. 3. This subproblem is recognised as mergeable and, therefore, we numerically merge the five clusters. 4. We merge both merged clusters CD6 to 10 and CD1 to 5 into a single cluster by employing Rule 3. Example. In the graph from Fig. 7c, the elimination process successively removes clusters CD5, CD8, CD7, CD6, CD10, and CD9. Then it determines equality NSE = 3(k  1) = 9 for the subproblem with k = 4 clusters CD1 to CD4. It does not find smaller mergeable subproblems. The further steps are: 1. We numerically merge the structurally well-constrained subproblem consisting of CD1 to CD4. 2. There are four inter-cluster constraints between the merged CD1 to 4 and CD5, as shown in Fig. 7d. We ignore one of them, and merge the clusters by employing Rule 3. Then we check the previously ignored structurally redundant constraint for consistency.

3. The merged CD1 to 4 (or eventually CD1 to 5) may be ignored during the analysis of the simplified graph since its valence is 3. We therefore deal with the subproblem consisting of CD6 to CD10. 4. We determine that this subproblem is structurally under-constrained. We numerically merge the six clusters CD1 to 4 (or eventually CD1 to 5) and CD6 to CD10. Example. We consider the problem from Fig. 4d. There are k = 8 clusters and 24 inter-cluster distances. The problem is supposed to be structurally over-constrained since NSE > 3(k  1) = 21. Each cluster has six neighbours and the valence 6. We eliminate CD1 and obtain a presumably structurally well-constrained subproblem with k = 7 clusters and NSE = 18 = 3(k  1). The further analysis does not find any mergeable subgraphs and, therefore, we numerically merge CD2 to CD8. A huge system of 24 equations has to be solved. After this, we have two clusters: CD1 and the merged CD2 to 8 connected by six inter-cluster constraints. Three of them are ignored and the problem is reduced by employing Rule 3. After this, we check the previously ignored redundant constraints for consistency.

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

5.4. Management of local coordinate frames Corollary 2 of Definition 4 states that a cluster is easily constructible by local propagation in linear time. This is certainly true for the initial clusters where, owing to preprocessing by calculation of redundant dimensions, all internal distances and angles are known. Each initial cluster CDj has its local polar coordinate frame with the origin in a selected cluster’s point called the anchor point anj, and with the polar axis towards another selected point called the azimuth point azj. Of course, the polar axis has no relevance in 1-clusters. The slope angle of the polar axis in the global (or local parent’s) coordinate frame will be referred to as the azimuth uj of CDj. Algorithm 6 represents a simple linear-time local propagation algorithm for positioning points of the initial cluster CD. It is also applicable for any other point set P = {pk}, 0 < k < |P|, with known polar points’ coordinates (rk, ak) in the attached local coordinate frame described by its anchor point an and its azimuth point az. Here rk represents the known distance | an pk| from the matrix of distances, and ak is the known angle \(az  an, pk  an) from the matrix of angles. Algorithm 6. Local Propagation(CD, an, az) Azimuth u :¼ slope angle of the vector (az – an). for (each point pk 2 CDn{an, az}) do pk :¼ an + rk(cos (u + ak), sin(u + ak)). A particular cluster merging operation uses the temporary points’ coordinates of a previously selected reference cluster to determine the temporary positions of the remaining involved clusters’ points. Once we have n points positioned, it is then easy to measure all O(n2) distances and all O(n4) angles. However, we try to optimise the algorithm and to compute only those dimensions that are really necessary. For this purpose, each merged cluster CDi is organised into a hierarchical data structure. A selected initial cluster, merged into CDi, named the anchor initial cluster of CDi, and labelled AN(CDi), becomes father of the remaining involved initial clusters. Among these remaining initial clusters, we select the azimuth initial cluster AN(CDi). Just like the incorporated initial clusters, the merged cluster CDi also has its local coordinate frame. The anchor point of CDi, labelled an(CDi), coincides with the anchor point of AN(CDi), and the azimuth point of CDi, az(CDi), is equal to the anchor point of AZ(CDi). The slope angle of the polar axis of CDi in the global (or local parent’s) coordinate frame is referred to as the azimuth u(CDi). When we merge clusters, we only have to define the transformations between the coordinate frames of the merged cluster CDi and the incorporated initial clusters CDi,j. Only three transformation parameters have to be determined for each CDi,j: (1) the displacement di,j :¼ |AN(CDi)ani,j|, (2) the displacement angle hi,j :¼ \(AZ(CDi)  AN(CDi), ani,j  AN(CDi)), or 0

783

when ani,j = AN(CDi), and (3) the already mentioned azimuth ui,j :¼ \(ani,j  AN(CDi), azi,j  ani,j). Theorem 8 assures that all n(n  1)/2 angles between n vectors in a CA set are known but, typically, the majority of them are never used. An optimisation can save a significant amount of time. The key of optimisation is the data structure. It actually only upgrades the hierarchical organisation of merged and initial clusters by an additional top level representing CA sets. Just like a cluster, a CA set CAi also has its anchor point an(CAi) and the azimuth point az(CAi) which coincide with the anchor point and the azimuth point of one of the included clusters. These points’ coordinates and the slope of the polar axis u(CAi), if not explicitly given, are read from the sketch. When we merge CA sets, only the azimuths u(CDj) of those clusters being moved into the extended CAi set must be updated. For each cluster CDj, we have to determine the displacement angle hi,(CDj), and then compute u(CDj) :¼ u(CAi) + hi,(CDj). Fig. 8 illustrates the meaning of the parameters of local coordinate frames of CA sets, merged and initial clusters. 5.5. Geometrical solvability Here we give some guidelines for dealing with geometrical redundancy. We consider this problem only in its most common form i.e. when too many determined angles appear in a loop of vectors. Theorem 4 partially solves this problem by allowing that, at most, a single angle links a pair of clusters in a particular cluster merging operation. This limitation must be dealt with when searching for minimal mergeable subproblems, or for patterns mergeable by any of Rules 3, 4, 6, 7 or 8. Theorem 4 prevents redundant angles between independently considered pairs of clusters. Theorem 13, on the other hand, prevents the use of too many angles between the clusters in the group of more than two clusters. Theorem 13. Let N 0 represent a set of clusters, and E 0 a set of angles between pairs of clusters from N 0 . The angles from E 0 are not independent if the graph G 0 = (N 0 , E 0 ) contains cycles. The proof is trivial. A loop in the graph G 0 corresponds to a polygon where all the internal angles are explicitly defined. These angles are certainly not independent. Theorem 13 must be taken into consideration before employing Rule 8 and during the heuristic search for minimal mergeable subproblems. It actually states that a cluster merging operation on m clusters may use, at most, m  1 angles between the clusters. More precisely, any subproblem with k 6 m clusters contains, at most, k  1 independent angles between the clusters. If the clusters involved in the merging operation belong to different CA sets, then we must separately fulfil this condition for each CA set. With other words, in a group of k clusters from the same CA set, angles between pairs of clusters contribute k  1 to the number of scalar equations NSE as already stated in Section 5.3.

784

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

Fig. 8. The hierarchical organisation of line sets, merged and initial clusters.

6. Results Here we demonstrate the use of the algorithm in several additional practical examples. The stress is on employing the construction rules, as some expressive examples of employing the heuristic search for the mergeable subgraph have been given in the previous section. Example. We consider the quadrilateral from Fig. 3a. Point p5 is inserted while transforming angles into the equations. We obtain three clusters as shown in Fig. 9a. Rule 2 merges CD2 and CD3, since they share point p3 and they are linked by e2  Angle(l2, l3, p  c). After this, Rule 3 merges CD1 and CD2,3, linked by three constraints e1  Distance(p1, p2, a), e3  Angle(l4, l3, d), and e4  Collinear(p3, p5, p4). Example. We consider the quadrilateral from Fig. 3b. Points p5 and p6 are added to enable transformation of

the angle between the diagonals into the equations. We obtain k = 4 clusters as shown in Fig. 9b. None of Rules 1–8 can be used and, therefore, we check the entire subgraph for structural rigidity. It is structurally well-constrained since NSE = 9 = 3(k  1) and because the algorithm does not find smaller mergeable subgraphs. The shared p1 and p4 each contribute 2 to NSE, and there are five inter-cluster constraints e1  Distance(p1, p2, a), e2  Distance(p3, p4, c), e3 defining the angle between the diagonals i.e. u(p1, p5, p2, p6)  cos u = 0, e4  Collinear(p2, p6, p4), and e5  Collinear(p1, p5, p3). The entire GCP is solved numerically. The frontier points p1 and p4 from CD1 are initially positioned. We need 8 equations to determine coordinates of the remaining four points. Besides e1 to e5, we use the internal cluster distances e6  Distance(p1, p5, 1), e7  Distance(p2, p3, b), and e8  Distance(p2, p6, 1).

Fig. 9. Clusters and inter-cluster constraints in the quadrilaterals from (a) Fig. 3a, and (b) from Fig. 3b.

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

785

Fig. 10. Quadrilateral constrained by three dimensional and two conditional constraints.

Example. We consider the quadrilateral p1p2p3p4 from Fig. 10a. Five constraints are given: Angle(l6, l3, d), Distance(p1, p3, e), Distance(p2, p4, f), Sum6(p1, p2, p2, p3, p3, p4), and Sum4L(l5, l6, l5, l1, a). The last one is illustrated by the angles u and u  a. The constraints Angle and Sum4L require the introduction of new points p5 and p6. The decomposition results in four clusters. Rule 2 immediately merges CD3 and CD4, since they share point p4 and they are linked by constraint Angle(l6, l3, d), as shown in Fig. 10b. None of Rules 1–8 can then reduce k = 3 clusters CD1, CD2, and CD3, 4. The GCP is structurally well-constrained since NSE = 6 = 3(k  1) and because the algorithm does not find smaller mergeable components. The shared p1 contributes 2 to NSE, and there are also four inter-cluster constraints: e1  Collinear(p1, p5, p2), e2  Collinear(p3, p4, p6), e3  Sum6(p1, p2, p2, p3, p3, p4), and e4  Sum4L(l5, l6, l5, l1, a). Note that these last two are represented by hyperedges, since each of them links three clusters. We solve this subproblem numerically. Local propagation positions points p2, p4 and p6 from CD3,4. We then need six equations to determine six coordinates of p1, p3, and p5. Besides the equations e1 to e4, we use the internal cluster distances e5  Distance(p1, p5, 1), and e6  Distance(p1, p3, e).

Example. The example in Fig. 11 demonstrates the use of CA sets. The dashed lines represent angles between the clusters. The initial problem in Fig. 11a is structurally under-constrained. In addition, none of Rules 1–8 can be directly employed. However, Theorem 10 allows transformation of Sum6L into the equivalent constraint Angle. This enables us to merge CA1 and CA2. Since all the angles in a CA set are known, we can determine the redundant angle between CD1 and CD3 (Fig. 11b). This angle enables the use of Rule 3 to merge CD1 and CD3, as shown in Fig. 11c. We then ignore one of the angles between CD1,3 and CD2 (use of Theorem 4), and merge the structurally under-constrained subproblem by employing Rule 8.

7. Summary This paper presents our original shape-recognitionbased constructive geometric constraint solver. It incorporates pre-processing by calculation of redundant dimensions, originally presented in [25], cluster creation, the original step of cluster registration, organisation of clusters into generally larger CA sets, the recombination phase with numerous original functionalities and, finally, an optimised

Fig. 11. The use of CA sets can enable us to avoid numerical solving of systems of equations.

786

D. Podgorelec et al. / Advances in Engineering Software 39 (2008) 770–786

local-propagation-based construction. The main original contributions presented in the paper are: transformation of a rich repertoire of constraints into the equations on points’ coordinates only (Section 3.1), the cluster registration step, assuring that the number of clusters for further consideration does not exceed the number of points n (Section 4.1), some special rules (Theorems 9–12) for production of redundant angles from angular conditional constraints, the heuristic approach for determining minimal mergeable subproblems (Section 5.3), the hierarchical organisation of CA sets, merged and initial clusters, which enables calculation of only those redundant dimensions that are really needed (Section 5.4), and consideration of geometrical redundancy caused by too many angles in a loop of vectors (Section 5.5). However, there are still numerous problems left for future investigation. The heuristic structural analysis requires further practical tests, theoretical confirmation or eventual correction of the selected heuristics. The preprocessing and the decomposition phase must be optimised. Finally, an incremental version of the presented algorithm is perhaps the most challenging. Acknowledgement The authors are grateful to the Slovenian Research Agency for supporting this research under the project P20041 Computer systems, methodologies and intelligent services. References [1] Hoffmann CM, Lomonosov A, Sitharam M. Geometric constraint decomposition. In: Bruderlin B, Roller D, editors. Geometric constraint solving and applications; 1998, p. 170–95. [2] Lamure H, Michelucci D. Solving geometric constraints by homotopy. IEEE Trans Visualization Comp Graphics 1996;2:28–34. [3] Ge J-X, Chou S-C, Gao X-S. Geometric constraint satisfaction using optimization methods. Comput-Aided Des 1999;31:867–79. [4] Kapur D, Mundy JL. Wu’s method and its application to perspective viewing. Artif Intell 1988;37:15–36. [5] Gao X-S, Chou S-C. Solving geometric constraint systems. II.A symbolic approach and decision of Rc-constructibility. ComputAided Des 1998;30:115–22. [6] Sridhar N, Agrawal R, Kinzel GL. Algorithms for the structural diagnosis and decomposition of sparse, underconstrained design systems. Comput-Aided Des 1996;28:237–49.

ˇ alik B, Guid N, Clapworthy G. Constraint-based object modelling. [7] Z J Eng Des 1996;7:209–32. [8] Owen J. Algebraic solution for geometry from dimensional constraints. In: ACM symposium on foundations of solid modeling, Austin (TX); 1991. p. 397–407. [9] Fudos I, Hoffmann CM. A graph-constructive approach to solving systems of geometric constraints. ACM Trans Graphics 1997;16: 179–216. [10] Joan-Arinyo R, Soto-Riera A, Vila-Marta S, Vilaplana-Pasto´ J. Revisiting decomposition analysis of geometric constraint graphs. ACM solid modeling 02. New York: ACM Press; 2003. [11] Lee K-Y, Kwon O-H, Lee J-Y, Kim T-W. A hybrid approach to geometric constraint solving with graph analysis and reduction. Adv Eng Soft 2003;34:103–13. [12] Hoffmann CM, Lomonosov A, Sitharam M. Decomposition plans for geometric constraint systems, Part I: performance measures for CAD. J Symb Comp 2001;31:367–408. [13] Hoffmann CM, Lomonosov A, Sitharam M. In: Smolka G, editor. Finding solvable subsets of constraint graphs. LNCS 1330. New York: Springer; 1997. p. 463–77. [14] Hoffmann CM, Lomonosov A, Sitharam M. Decomposition plans for geometric constraint systems, Part II: new algorithms. J Symb Comp 2001;31:409–27. [15] Gao X-S, Zhang G-F. Geometric constraint solving via c-tree decomposition. ACM solid modeling 02. New York: ACM Press; 2003, p. 45–55. [16] Latham RS, Middletich AE. Connectivity analysis: a tool for processing geometric constraints. Comput-Aided Des 1996;28: 917–28. [17] Sitharam M, Zhou Y. A tractable, approximate, combinatorial 3d rigidity characterization. In: Proceedings of automated deduction in geometry. Gainesville, USA: University of Florida; 2004. [18] Foufou S, Michelucci D, Jurzak J-P. Numerical decomposition of geometric constraints. In: ACM symposium on solid and physical modelling. Boston, USA: MIT, ACM Press; 2005. p. 143–51. [19] Michelucci D, Foufou S. Geometric constraint solving: the witness configuration method. Comput-Aided Des 2006;38:284–99. [20] Jermann C, Neveu B, Trombettoni G. Algorithms for identifying rigid subsystems in geometric constraint systems. In: IJCAI’03, Acapulco, Mexico; 2003. p. 233–8. [21] Pabon JA. Modeling method for sorting dependencies among geometric entities. US States Patent 5, 251, 290; October 1993. [22] Sunde G. A CAD system with declarative specification of shape. Eurographics workshop on intelligent CAD systems, Noordwijkerhout, The Netherlands; 1987. p. 90–104. [23] Verroust A, Schonek F, Roller D. Rule-oriented method for parameterized computer-aided design. Comput-Aided Des 1992;24: 209–14. [24] Hoffmann CM, Sitharam M, Yuan B. Making constraint solvers more usable: overconstraint problem. Comput-Aided Des 2004;36: 377–99. [25] Podgorelec D. A new constructive approach to constraint-based geometric design. Comput-Aided Des 2002;34:769–857.