added simple examples for the current cpp implementations

Author: tomvanmele

src/compas/geometry/algorithms/_smoothing_cpp/smoothing_cpp_test.py

@@ -0,0 +1,70 @@
+import compas
+from compas.datastructures import Mesh
+from compas.plotters import MeshPlotter
+from compas.geometry import smooth_centroid_cpp
+
+kmax = 50
+
+# make a mesh
+# and set the default vertex and edge attributes
+
+mesh = Mesh.from_obj(compas.get('faces.obj'))
+
+# extract numerical data from the datastructure
+
+vertices  = mesh.get_vertices_attributes(('x', 'y', 'z'))
+adjacency = [mesh.vertex_neighbours(key) for key in mesh.vertices()]
+fixed     = [int(mesh.vertex_degree(key) == 2) for key in mesh.vertices()]
+
+slider = list(mesh.vertices_where({'x': (-0.1, 0.1), 'y': (9.9, 10.1)}))[0]
+
+# make a plotter for (dynamic) visualization
+# and define a callback function
+# for plotting the intermediate configurations
+
+plotter = MeshPlotter(mesh, figsize=(10, 7))
+
+def callback(k, xyz):
+    print(k)
+
+    if k < kmax - 1:
+        xyz[slider][0] = 0.1 * (k + 1)
+
+    plotter.update_vertices()
+    plotter.update_edges()
+    plotter.update(pause=0.001)
+
+    for key, attr in mesh.vertices(True):
+        attr['x'] = xyz[key][0]
+        attr['y'] = xyz[key][1]
+        attr['z'] = xyz[key][2]
+
+# plot the lines of the original configuration of the mesh
+# as a reference
+
+lines = []
+for u, v in mesh.edges():
+    lines.append({
+        'start': mesh.vertex_coordinates(u, 'xy'),
+        'end'  : mesh.vertex_coordinates(v, 'xy'),
+        'color': '#cccccc',
+        'width': 0.5
+    })
+
+plotter.draw_lines(lines)
+
+# draw the vertices and edges in the starting configuration
+# and pause for a second before starting the dynamic visualization
+
+plotter.draw_vertices(facecolor={key: '#000000' for key in mesh.vertices() if mesh.vertex_degree(key) == 2})
+plotter.draw_edges()
+
+plotter.update(pause=0.5)
+
+# run the smoother
+
+xyz = smooth_centroid_cpp(vertices, adjacency, fixed, kmax=kmax, callback=callback)
+
+# keep the plot alive
+
+plotter.show()

src/compas/geometry/algorithms/_smoothing_cpp/smoothing_cpp_test_rhino.py

@@ -0,0 +1,78 @@
+from __future__ import print_function
+from __future__ import absolute_import
+from __future__ import division
+
+from ctypes import c_double
+import compas
+import compas_rhino
+
+from compas.datastructures import Mesh
+from compas.geometry import smooth_centroid_cpp
+
+from compas_rhino.helpers import MeshArtist
+from compas_rhino.conduits import PointsConduit
+from compas_rhino.conduits import LinesConduit
+
+kmax = 50
+
+# make a mesh
+# and set the default vertex and edge attributes
+
+mesh = Mesh.from_obj(compas.get('faces.obj'))
+
+edges = list(mesh.edges())
+
+# extract numerical data from the datastructure
+
+vertices  = mesh.get_vertices_attributes(('x', 'y', 'z'))
+adjacency = [mesh.vertex_neighbours(key) for key in mesh.vertices()]
+fixed     = [int(mesh.vertex_degree(key) == 2) for key in mesh.vertices()]
+
+# make a artist for (dynamic) visualization
+# and define a callback function
+# for plotting the intermediate configurations
+
+slider = list(mesh.vertices_where({'x': (-0.1, 0.1), 'y': (9.9, 10.1)}))[0]
+
+artist = MeshArtist(mesh, layer='SmoothMesh')
+artist.clear_layer()
+
+
+def callback(k, xyz):
+    compas_rhino.wait()
+
+    print(k)
+
+    if k < kmax - 1:
+        xyz[slider][0] = c_double(0.1 * (k + 1))
+
+    pointsconduit.points = [mesh.vertex_coordinates(key) for key in mesh.vertices()]
+    linesconduit.lines = [mesh.edge_coordinates(u, v) for u, v in edges]
+
+    pointsconduit.redraw(k)
+    linesconduit.redraw(k)
+
+    for key, attr in mesh.vertices(True):
+        attr['x'] = float(xyz[key][0])
+        attr['y'] = float(xyz[key][1])
+        attr['z'] = float(xyz[key][2])
+
+
+pointsconduit = PointsConduit(radius=10, refreshrate=5)
+linesconduit = LinesConduit(refreshrate=5)
+
+
+with pointsconduit.enabled():
+    with linesconduit.enabled():
+        xyz = smooth_centroid_cpp(vertices, adjacency, fixed, kmax=kmax, callback=callback)
+
+
+for key, attr in mesh.vertices(True):
+    attr['x'] = xyz[key][0]
+    attr['y'] = xyz[key][1]
+    attr['z'] = xyz[key][2]
+
+artist.clear_edges()
+artist.draw_vertices()
+artist.draw_edges()
+artist.redraw()

src/compas/geometry/algorithms/_smoothing_cpp/src/main.cpp

@@ -3,8 +3,6 @@
 // cd _smoothing_cpp
 // g++ -shared -fPIC src/main.cpp -o smoothing.so
 
-using namespace std;
-
 extern "C"
 {
     typedef void callback(int k);
@@ -19,7 +17,7 @@ void smooth_centroid(int v, int *nbrs, int *fixed, double **vertices, int **neig
     int j, n;
     double cx, cy, cz;
 
-    vector< vector<double> > xyz(v, vector<double>(3, 0.0));
+    std::vector< std::vector<double> > xyz(v, std::vector<double>(3, 0.0));
 
     for (k = 0; k < kmax; k++) {
 

src/compas/numerical/algorithms/_fd_cpp/fd.dll

@@ -0,0 +1,44 @@
+from __future__ import print_function
+from __future__ import absolute_import
+from __future__ import division
+
+from ctypes import c_double
+import random
+import compas
+import compas_rhino
+
+from compas.datastructures import Network
+from compas.numerical import fd_cpp
+from compas.utilities import i_to_red
+
+from compas_rhino.helpers import NetworkArtist
+
+network = Network.from_obj(compas.get('saddle.obj'))
+
+# for u, v, attr in network.edges(True):
+#     attr['q'] = 1.0 * random.randint(1, 5)
+
+vertices = network.get_vertices_attributes('xyz')
+loads    = network.get_vertices_attributes(('px', 'py', 'pz'), (0.0, 0.0, 0.0))
+fixed    = network.leaves()
+free     = [i for i in range(len(vertices)) if i not in fixed]
+edges    = list(network.edges())
+q        = network.get_edges_attribute('q', 1.0)
+
+
+artist = NetworkArtist(network, layer='Network')
+
+artist.clear_layer()
+artist.draw_edges(color='#cccccc')
+artist.redraw()
+
+xyz = fd_cpp(vertices, edges, fixed, q, loads)
+
+for key, attr in network.vertices(True):
+    attr['x'] = xyz[key][0]
+    attr['y'] = xyz[key][1]
+    attr['z'] = xyz[key][2]
+
+artist.draw_vertices(color={key: '#ff0000' for key in network.leaves()})
+artist.draw_edges(color='#000000')
+artist.redraw()

src/compas/numerical/algorithms/_fd_cpp/fd.so

@@ -0,0 +1,91 @@
+from __future__ import print_function
+from __future__ import absolute_import
+from __future__ import division
+
+import Rhino
+
+from System.Drawing.Color import FromArgb
+
+from Rhino.Geometry import Point3d
+from Rhino.Geometry import Line
+
+import compas
+import compas_rhino
+
+from compas.datastructures import Network
+from compas.numerical import fd_cpp
+
+from compas_rhino.helpers import NetworkArtist
+from compas_rhino.helpers.selectors import VertexSelector
+from compas_rhino.conduits import LinesConduit
+
+
+network = Network.from_obj(compas.get('saddle.obj'))
+
+
+vertices = network.get_vertices_attributes('xyz')
+loads    = network.get_vertices_attributes(('px', 'py', 'pz'), (0.0, 0.0, 0.0))
+fixed    = network.leaves()
+free     = [i for i in range(len(vertices)) if i not in fixed]
+edges    = list(network.edges())
+q        = network.get_edges_attribute('q', 1.0)
+
+
+artist = NetworkArtist(network, layer='Network')
+
+artist.clear_layer()
+artist.draw_edges(color='#cccccc')
+artist.redraw()
+
+xyz = fd_cpp(vertices, edges, fixed, q, loads)
+
+for key, attr in network.vertices(True):
+    attr['x'] = xyz[key][0]
+    attr['y'] = xyz[key][1]
+    attr['z'] = xyz[key][2]
+
+artist.draw_vertices(color={key: '#ff0000' for key in network.leaves()})
+artist.draw_edges(color='#000000')
+artist.redraw()
+
+# select and drag one of the anchors
+# update equilibrium using *OnDynamicDraw*
+
+move = VertexSelector.select_vertex(network)
+move = int(move)
+
+sp = network.vertex_coordinates(move)
+ep = sp[:]
+ep[2] += 1
+
+vertical = Line(Point3d(*sp), Point3d(*ep))
+
+color = FromArgb(255, 255, 255)
+
+
+def OnDynamicDraw(sender, e):
+    global xyz
+    xyz[move] = list(e.CurrentPoint)
+    xyz = fd_cpp(xyz, edges, fixed, q, loads)
+    for u, v in edges:
+        e.Display.DrawLine(Point3d(* xyz[u]), Point3d(* xyz[v]), color)
+
+
+gp = Rhino.Input.Custom.GetPoint()
+gp.SetCommandPrompt('Point to move to?')
+gp.Constrain(vertical)
+gp.DynamicDraw += OnDynamicDraw
+gp.Get()
+
+if gp.CommandResult() == Rhino.Commands.Result.Success:
+    xyz[move] = list(gp.Point())
+    xyz = fd_cpp(xyz, edges, fixed, q, loads)
+
+    for key, attr in network.vertices(True):
+        attr['x'] = xyz[key][0]
+        attr['y'] = xyz[key][1]
+        attr['z'] = xyz[key][2]
+
+    artist.draw_vertices(color={key: '#ff0000' for key in network.leaves()})
+    artist.draw_edges(color='#000000')
+    artist.redraw()

src/compas/numerical/algorithms/_fd_cpp/fd_cpp_test_rhino.py

@@ -6,10 +6,10 @@ extern "C"
 {
 	typedef void callback(int k);
 
-	void fd(int numv, int nume, int numfix, double **vertices, int **edges, double **loads, double *q, int *fixed, int *free);
+	void fd(int numv, int nume, int numfix, double** vertices, int** edges, double** loads, double* q, int* fixed, int* free);
 }
 
-void fd(int numv, int nume, int numfix, double **vertices, int **edges, double **loads, double *q, int *fixed, int *free)
+void fd(int numv, int nume, int numfix, double** vertices, int** edges, double** loads, double* q, int* fixed, int* free)
 {
 	int i;
 	int numfree = numv - numfix;
@@ -70,8 +70,8 @@ void fd(int numv, int nume, int numfix, double **vertices, int **edges, double *
 
 	Cit = Ci.transpose();
 
-	A = Ci.transpose() * Q * Ci;
-	b = Pi - Cit * Q * Cf * Xf;
+	A.noalias() = Cit * Q * Ci;
+	b.noalias() = Pi - Cit * Q * Cf * Xf;
 
 	Xi = A.colPivHouseholderQr().solve(b);
 

src/compas/numerical/algorithms/_fd_cpp/fd_cpp_test_rhino_interactive.py

@@ -37,12 +37,14 @@ def fd_cpp(vertices, edges, fixed, q, loads):
             print(DLL)
             raise
 
+    free = [i for i in range(len(vertices)) if i not in fixed]
+
     cvertices = Array2D(vertices, 'double')
-    cedges = Array2D(edges, 'int')
-    cloads = Array2D(loads, 'double')
-    cq = Array1D(q, 'double')
-    cfixed = Array1D(fixed, 'int')
-    cfree = Array1D(free, 'int')
+    cedges    = Array2D(edges, 'int')
+    cloads    = Array2D(loads, 'double')
+    cq        = Array1D(q, 'double')
+    cfixed    = Array1D(fixed, 'int')
+    cfree     = Array1D(free, 'int')
 
     lib.fd.argtypes = [
         ctypes.c_int,
@@ -86,7 +88,6 @@ def fd_cpp(vertices, edges, fixed, q, loads):
     vertices = network.get_vertices_attributes('xyz')
     loads    = network.get_vertices_attributes(('px', 'py', 'pz'), (0.0, 0.0, 0.0))
     fixed    = network.leaves()
-    free     = [i for i in range(len(vertices)) if i not in fixed]
     edges    = list(network.edges())
     q        = network.get_edges_attribute('q', 1.0)