added mapping for lines

TODO.txt

@@ -3,7 +3,7 @@ TODO/WISH LIST
 
 Code: Some bits of the code are ugly:
     - vector module, matpy module and Numpy package, all offer similar
-      functionality. Use only one, get rid of two. 
+      functionality. Use only one (numpy), get rid of two. 
     - more documentation!
 
 Speed: The solver is currently much too slow. The problem is the pattern

geosolver/clsolver.py

@@ -590,7 +590,7 @@ class ClusterSolver(Notifier):
             if isinstance(item, Cluster):
                 for var in item.vars:
                     if len(self.find_dependend(var)) == 0:
-                        self._graph.remove(var)
+                        self._graph.rem_vertex(var)
             # notify listeners
             self.send_notify(("remove", item))
         # restore toplevel (also added to _new)

geosolver/configuration.py

@@ -333,6 +333,11 @@ class Configuration:
         else:
             return "Configuration("+str(self.map)+")"
 
+    def __contains__(self,var):
+        return var in self.map
+
+    def __getitem__(self,var):
+        return self.map[var]
 
 def testeq():
     p1 = vector.vector([0.0,0.0,0.0])

geosolver/geometric.py

@@ -7,13 +7,14 @@ from clsolver import PrototypeMethod, SelectionMethod
 from clsolver3D import ClusterSolver3D 
 from clsolver2D import ClusterSolver2D 
 from cluster import *
+from selconstr import SelectionConstraint
 from configuration import Configuration 
 from diagnostic import diag_print
 from constraint import Constraint, ConstraintGraph
 from notify import Notifier, Listener
 from tolerance import tol_eq
 from intersections import angle_3p, distance_2p
-from selconstr import SelectionConstraint
+from intersections import distance_point_line
 from intersections import is_left_handed, is_right_handed
 from intersections import is_clockwise, is_counterclockwise 
 from intersections import transform_point, make_hcs_3d
@@ -107,7 +108,7 @@ class GeometricProblem (Notifier, Listener):
         if variable in self.prototype:
             return self.prototype[variable]
         else:
-            raise StandardError, "unknown variable variable"
+            raise StandardError, "unknown variable "+str(variable)
 
     # --------------- points - depricated  -------------
     def add_point(self, variable, prototype):
@@ -201,18 +202,39 @@ class GeometricProblem (Notifier, Listener):
         else:
             return None
 
+    def get_constraints_with_type_on_variables(self, constrainttype, variables):
+        candidates = None
+        for var in variables:
+            if candidates == None:
+                candidates = set(filter(lambda c:isinstance(c,constrainttype),self.cg.get_constraints_on(var)))
+            else:    
+                candidates.intersection_update(filter(lambda c:isinstance(c,constrainttype),self.cg.get_constraints_on(var)))
+        return candidates
+
+    def get_unique_constraint(self, constrainttype, variables):
+        candidates = self.get_constraints_with_type_on_variables(constrainttype, variables)
+        if len(candidates) == 0:
+            return None
+        elif len(candidates) == 1:
+            return candidates[0]
+        else: # >= 1
+            raise StandardError, "multiple constraints found"
+        
     def get_coincidence(self, p, g):
-        print "GeometricProblem.get_coincidence NOT IMPLEMENTED"
+        return self.get_unique_constraint(CoincidenceConstraint, [p,g])
-        return None
 
-    def get_rigid(self, vars):
+    def get_rigid(self, variables):
-        print "GeometricProblem.get_rigid NOT IMPLEMENTED"
+        return self.get_unique_constraint(RigidConstraint, variables)
-        return None
 
-    def get_mate(self, vars):
+    def get_mate(self, variables):
-        print "GeometricProblem.get_mate NOT IMPLEMENTED"
+        return self.get_unique_constraint(MateConstraint, variables)
-        return None
 
+    def get_coincident_points(self, geometry):
+        coincidences = self.get_constraints_with_type_on_variables(CoincidenceConstraint, [geometry])
+        points = set()
+        for constraint in coincidences:
+            points.update(filter(lambda var: isinstance(var, Point) and var != geometry, constraint.variables()))
+        return points
 
     def verify(self, solution):
         """returns true iff all constraints satisfied by given solution. 
@@ -384,12 +406,13 @@ class GeometricSolver (Listener):
             drclusters = map[geocluster]
             drcluster = max(drclusters, key=lambda c: c.creationtime)
             # determine solutions
-            solutions = self.dr.get(drcluster)
+            geocluster.solutions = self._map_cluster_solutions(drcluster)
+            # determine incidental underconstrainedness
             underconstrained = False
-            if solutions != None:
+            configurations = self.dr.get(drcluster)
-                for solution in solutions:
+            if configurations != None:
-                    geocluster.solutions.append(solution.map)
+                for config in configurations:
-                    if solution.underconstrained:
+                    if config.underconstrained:
                         underconstrained = True
             # determine flag
             if drcluster.overconstrained:
@@ -442,14 +465,51 @@ class GeometricSolver (Listener):
     def get_solutions(self):
         """Returns a list of Configurations, which will be empty if the
            problem has no solutions. Note: this method is
-           cheaper but less informative than get_decomposition. The
+           cheaper but less informative than get_decomposition. 
-           list and the configurations should not be changed (since they are
+        """
-           references to objects in the solver)."""
+        #"""The list and the configurations should not be changed (since they are
+        #references to objects in the solver)."""
+        # find top level rigid and all its configurations
         rigids = filter(lambda c: isinstance(c, Rigid), self.dr.top_level())
         if len(rigids) != 0:   
-            return self.dr.get(rigids[0])
+            solutions = self._map_cluster_solutions(drcluster)
         else:
-            return []
+            solutions = []
+        return solutions
+
+    def _map_cluster_solutions(self, drcluster):
+        # map dr-cluster configurations to solutions, i.e. a map from problem variables to values           
+        configurations = self.dr.get(drcluster)
+        solutions = []
+        diag_print("mapping cluster "+str(drcluster)+" #configurations="+str(len(configurations)),"GeometricSolver")
+        for configuration in configurations:
+            solution = {}
+            for var in self.problem.cg.variables():
+                if isinstance(var, Point):
+                    assert len(self._map[var].vars) == 1
+                    point = iter(self._map[var].vars).next()
+                    if point in configuration:
+                        solution[var] = configuration[point]
+                elif isinstance(var, Line):
+                    if var in self._map:
+                        vertices = list(self._map[var].vars) 
+                    else: 
+                        # when line coincident with 2 points, then not mapped to cluster... use any two coincident points
+                        points = self.problem.get_coincident_points(var)
+                        assert len(points) >= 2
+                        print "points",points
+                        vertices = map(lambda v: iter(self._map[v].vars).next(), list(points)[0:2])
+                        print "vertices",vertices
+                    assert len(vertices) == 2
+                    if vertices[0] in configuration and vertices[1] in configuration:
+                        solution[var] = vector.vector(configuration[vertices[0]]).concatonated( vector.vector(configuration[vertices[1]]) )
+                        
+                else:
+                    raise StandardError, "unknown variable type"
+            #for
+            solutions.append(solution)
+        #for    
+        return solutions
 
     def get_status(self):
         """Returns a symbolic flag, one of:
@@ -543,6 +603,49 @@ class GeometricSolver (Listener):
             self.dr.add(rigid)
             self._update_variable(var)
     
+    def _add_line(self, var):
+        diag_print("add line "+str(var),"GeometricSolver")
+        # find coincident points
+        points = list(self.problem.get_coincident_points(var))
+        diag_print("on "+str(points),"GeometricSolver")
+        if len(points) == 0:
+            self._map_line_distance(var)
+        elif len(points) == 1:
+            self._map_line_point_distance(var, points[0])
+        elif len(points) == 2:
+            self._map_line_point_point(var, points[0], points[1])
+        else: # >=3
+            self._map_line_points_radial(var, points)
+    
+    def _map_line_distance(self,line):
+        # map a line (coincident with no points) to a distance cluster (on two new point variables)
+        v1 = str(line)+"_v1"
+        v2 = str(line)+"_v2"
+        dist = Rigid([v1,v2])
+        self._map[line] = dist
+        self._map[dist] = line
+        self.dr.add(dist)
+        self._update_variable(line)
+ 
+    def _map_line_point_distance(self,line, point):
+        # map a line coincident with one point to a distance clusters (and one new point variable)
+        v1 = list(self._map[point].vars)[0]
+        v2 = str(line)+"_v2"
+        dist = Rigid([v1,v2])
+        self._map[line] = dist
+        self._map[dist] = line
+        self.dr.add(dist)
+        self._update_variable(line)
+ 
+    def _map_line_point_point(self, line, point1, point2):
+        # map a line coincident with two existing point variables; no new clusters created
+        self._update_variable(line)
+
+    def _map_line_points_radial(self, line, points):
+        # map a line coincient with three or more points to a radial clusters 
+        raise NotImplementedError
+
+ 
     def _rem_variable(self, var):
         diag_print("GeometricSolver._rem_variable","gcs")
         if var in self._map:
@@ -600,6 +703,14 @@ class GeometricSolver (Listener):
         elif isinstance(con, SelectionConstraint):
             # add directly to clustersolver
             self.dr.add_selection_constraint(con)
+        elif isinstance(con, CoincidenceConstraint):
+            # re-map lines, planes, etc
+            lines = filter(lambda var: isinstance(var,Line),con.variables())
+            if len(lines)==1:
+                line = iter(lines).next()
+                self._rem_variable(line)
+                self._add_line(line)
+            #endif
         else:
             raise StandardError, "unknown constraint type"
             pass
@@ -710,12 +821,41 @@ class GeometricSolver (Listener):
         else:
             raise StandardError, "unknown constraint type"
     
-    def _update_variable(self, variable):
+    def _update_variable(self, var):
+        if isinstance(var, Point):
+            self._update_point(var)
+        elif isinstance(var, Line):
+            self._update_line(var)
+        elif isinstance(var, Plane):
+            self._update_plane(var)
+        else:
+            # assume point - depricated
+            self._update_point(var)
+
+    def _update_point(self, variable):
         cluster = self._map[variable]
-        proto = self.problem.get_point(variable)
+        proto = self.problem.get_prototype(variable)
         conf = Configuration({variable:proto})
         self.dr.set(cluster, [conf])
 
+    def _update_line(self, variable):
+        # note: line may not be mapped to a cluster at all!
+        if variable in self._map:
+            cluster = self._map[variable]
+            proto = self.problem.get_prototype(variable)
+            vertices = list(cluster.vars);
+            v1 = vertices[0]
+            v2 = vertices[1]
+            assert self.problem.dimension==2 or self.problem.dimension==3
+            if self.problem.dimension == 2:
+                p1 = proto[0:2]
+                p2 = proto[2:4]
+            elif self.problem.dimension == 3:
+                p1 = proto[0:3]
+                p2 = proto[3:6]
+            conf = Configuration({v1:p1, v2:p2})
+            self.dr.set(cluster, [conf])
+
     def _update_fix(self):
         if self.fixcluster:
             vars = self.fixcluster.vars
@@ -826,7 +966,7 @@ class GeometricVariable:
         return hash(self.name)
    
     def __repr__(self):
-        return repr(self.__class__)+"("+repr(self.name)+")"
+        return self.__class__.__name__+"("+repr(self.name)+")"
     
 
 class Point(GeometricVariable):
@@ -1093,17 +1233,17 @@ class CoincidenceConstraint(Constraint):
 
     def satisfied(self, mapping):
         """return True iff mapping from variable names to points satisfies constraint""" 
-        if isinstance(geometry, Point):
+        if isinstance(self._geometry, Point):
             p1 = mapping[self._point]
             p2 = mapping[self._geometry]
             return tol_eq(distance_2p(p1,p2),0)   
-        elif isinstance(geometry, Line):
+        elif isinstance(self._geometry, Line):
             p = mapping[self._point]
             l = mapping[self._geometry]
             p1 = l[0:3]
             p2 = l[3:6]
             return tol_eq(distance_point_line(p, p1, p2),0)   
-        elif isinstance(geometry, Plane):
+        elif isinstance(self._geometry, Plane):
             p = mapping[self._point]
             l = mapping[self._geometry]
             p1 = l[0:3]

geosolver/intersections.py

@@ -231,6 +231,20 @@ def distance_2p(p1, p2):
     """
     return vector.norm(p2 - p1)
 
+def distance_point_line(p,l1,l2):
+    """distance from point p to line l1-l2"""
+    # v,w is p, l2 relative to l1
+    v = p-l1
+    w = l2-l1
+    # x = projection v on w
+    l = (vector.norm(v) * vector.norm(w))
+    if tol_eq(l,0):
+        x = 0*v
+    else:
+        x = v * vector.dot(v,w) / l
+    # result is distance x,v
+    return vector.norm(x-v)
+
 # ------ 2D
 
 def is_clockwise(p1,p2,p3):

geosolver/vector.py

@@ -98,6 +98,10 @@ class vector(list):
             """
             print self
 
+        def concatonated(self,other):
+            """this vector concatonated with another"""
+            return vector(list(self)+list(other))
+
 ###############################################################################
 
 

test/test_geometry.py

@@ -10,16 +10,35 @@ from geosolver.geometric import Point, Line, CoincidenceConstraint
 from geosolver.vector import vector 
 from geosolver.diagnostic import diag_select, diag_print
 
-def line_problem():
+def line_problem1():
     """A problem with a Point, a Line and a CoincicentConstraint"""
     problem = GeometricProblem(dimension=3)
-    problem.add_variable(Point('p1'),vector([0.0, 0.0, 0.0]))
+    problem.add_variable(Line('l1'),vector([0.0, 0.0, 0.0, 1.0, 1.0, 1.0]))
+    return problem 
+
+
+def line_problem2():
+    """A problem with a Point, a Line and a CoincicentConstraint"""
+    problem = GeometricProblem(dimension=3)
+    problem.add_variable(Point('p1'),vector([3.0, 2.0, 1.0]))
     problem.add_variable(Line('l1'),vector([0.0, 0.0, 0.0, 1.0, 1.0, 1.0]))
     problem.add_constraint(CoincidenceConstraint(Point('p1'), Line('l1')))
     return problem 
 
+
+def line_problem3():
+    """A problem with a Point, a Line and a CoincicentConstraint"""
+    problem = GeometricProblem(dimension=3)
+    problem.add_variable(Point('p1'),vector([3.0, 2.0, 1.0]))
+    problem.add_variable(Point('p2'),vector([1.0, 1.0, 1.0]))
+    problem.add_variable(Line('l1'),vector([0.0, 0.0, 0.0, 1.0, 1.0, 1.0]))
+    problem.add_constraint(CoincidenceConstraint(Point('p1'), Line('l1')))
+    problem.add_constraint(CoincidenceConstraint(Point('p2'), Line('l1')))
+    problem.add_constraint(DistanceConstraint(Point('p1'), Point('p2'), 5.0))
+    return problem 
+
 def test_line():
-    test(line_problem())
+    test(line_problem3())
 
 if __name__ == "__main__": 
     test_line()

workbench/decompositionView.py

@@ -138,6 +138,13 @@ class DecompositionView(QtGui.QDialog):
 
     def optimiseGraphLayout(self):
         print "optimising graph layout..."
+        # force due to overlapping overlaps 
+        force_cluster = 0.2
+        # force due to connection length
+        force_connection = 0.025
+        # force due to clusters overlapping connections
+        force_cluster_connection = 0.05
+
         # create a graph of clusters and connections
         graph = geosolver.graph.Graph()
         if self.ui.graphicsScene != None:
@@ -150,7 +157,7 @@ class DecompositionView(QtGui.QDialog):
        
         l = list(graph.vertices())
 
-        # iteratively implove layout 
+        # iteratively improve layout 
         for i in range(100):
             # clear forces 
             for c in l:
@@ -178,9 +185,11 @@ class DecompositionView(QtGui.QDialog):
                         norm =  numpy.linalg.norm(direction)
                         if norm != 0:
                             direction = direction / numpy.linalg.norm(direction)
+                        else:
+                            direction = numpy.array([0,0])
                         #direction[1] = 0.0
-                        c1.force += -force*direction / 10.0;
+                        c1.force += -force*direction * force_cluster;
-                        c2.force += force*direction / 10.0;
+                        c2.force += force*direction * force_cluster;
                         #print "force 1", c1.force
                         #print "force 2", c2.force
             
@@ -190,30 +199,43 @@ class DecompositionView(QtGui.QDialog):
                 c2 = e[1] 
                 box1 = c1.boundingRect().translated(c1.position)
                 box2 = c2.boundingRect().translated(c2.position)
+                # force 1: pull together on x
                 centerdiff = box2.center()-box1.center()
-                direction = numpy.array([centerdiff.x(),centerdiff.y()])
+                direction = numpy.array([centerdiff.x(),0])
                 norm =  numpy.linalg.norm(direction)
                 if norm != 0:
                     direction = direction / numpy.linalg.norm(direction)
-                #goal = max(box1.width(),box2.width())
+                else:
-                goal = 5 * box1.height()
+                    direction = numpy.array([0,0])
-                force = (norm - goal) / 50.0
+                goal = 0 
-                #direction[1] = 0.0
+                force = (norm - goal) * force_connection;
+                c1.force += +force*direction;
+                c2.force += -force*direction;
+                # force 2: keep y at distance and in layer order
+                direction = numpy.array([0, centerdiff.y()])
+                norm =  numpy.linalg.norm(direction)
+                if norm != 0:
+                    direction = direction / numpy.linalg.norm(direction)
+                else:
+                    direction = numpy.array([0,0])
+                goal = box1.height() + box2.height()
+                force = (norm - goal) * force_connection;
                 c1.force += +force*direction;
                 c2.force += -force*direction;
                 #print "force ", force
                 
+
             # determine forces due to clusters overlapping connections
             n = len(l)
             for c in graph.vertices():
                 for e in graph.edges():
                     box1 = c.boundingRect().translated(c.position)
-                    box1.setWidth(2*box1.width())
+                    box1.setWidth(0.5*box1.width())
-                    box1.setHeight(2*box1.height())
+                    box1.setHeight(0.5*box1.height())
                     con = graph.get(e[0],e[1])
                     box2 = con.boundingRect()
-                    box2.setWidth(2*box2.width())
+                    box2.setWidth(0.5*box2.width())
-                    box2.setHeight(2*box2.height())
+                    box2.setHeight(0.5*box2.height())
                     #print "box 1", box1 
                     #print "box 2", box2 
                     if box1.intersects(box2):
@@ -224,23 +246,26 @@ class DecompositionView(QtGui.QDialog):
                         norm =  numpy.linalg.norm(direction)
                         if norm != 0:
                             direction = direction / numpy.linalg.norm(direction)
-                        #direction[1] = 0.0
+                        else:
-                        c.force += -force*direction / 50.0;
+                            direction = numpy.array([0,0])
+                        c.force += -force*direction * force_cluster_connection;
+                        e[0].force += force*direction * force_cluster_connection;
+                        e[1].force += force*direction * force_cluster_connection;
                         #print "force 1", c1.force
                         #print "force 2", c2.force
 
-
             # apply forces 
             for c in l:
                 move = QtCore.QPointF(c.force[0],c.force[1])
                 c.position += move
                 c.translate(move.x(), move.y())
-        # done iterating
 
             # uppate connectors
             for e in graph.edges():
                 connector = graph.get(e[0],e[1])
                 connector.determinePath()
+
+        # done iterating
         print "done"
 
     def updateViewports(self):