Vectorizedect

.github/workflows/python-package.yml

@@ -16,7 +16,7 @@ jobs:
     strategy:
       fail-fast: false
       matrix:
-        python-version: ["3.9", "3.10", "3.11"]
+        python-version: ["3.10", "3.11"]
 
     steps:
     - uses: actions/checkout@v4

Extra_Notebooks/Mesh.ipynb

@@ -0,0 +1,79 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "1 0\n",
+      "4 5\n",
+      "9 8\n",
+      "10 1\n"
+     ]
+    }
+   ],
+   "source": [
+    "from CGAL.CGAL_Kernel import Point_2\n",
+    "from CGAL.CGAL_Triangulation_2 import Triangulation_2\n",
+    "from CGAL.CGAL_Triangulation_2 import Triangulation_2_Vertex_circulator\n",
+    "from CGAL.CGAL_Triangulation_2 import Triangulation_2_Vertex_handle\n",
+    "\n",
+    "points=[]\n",
+    "points.append( Point_2(1,0) )\n",
+    "points.append( Point_2(3,2) )\n",
+    "points.append( Point_2(4,5) )\n",
+    "points.append( Point_2(9,8) )\n",
+    "points.append( Point_2(7,4) )\n",
+    "points.append( Point_2(5,2) )\n",
+    "points.append( Point_2(6,3) )\n",
+    "points.append( Point_2(10,1) )\n",
+    "  \n",
+    "t=Triangulation_2()\n",
+    "t.insert(points)\n",
+    "\n",
+    "vc = t.incident_vertices(t.infinite_vertex())\n",
+    "\n",
+    "if vc.hasNext():\n",
+    "  done = vc.next();\n",
+    "  iter=Triangulation_2_Vertex_handle()\n",
+    "  while(1):\n",
+    "    iter=vc.next()\n",
+    "    print (iter.point())\n",
+    "    if iter == done:\n",
+    "      break\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "base",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.11.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}

Makefile

@@ -13,7 +13,7 @@ reqs:
 
 clean:
 	# Running autopep8
-	@autopep8 -r --in-place ect/
+	@autopep8 -r --in-place src/ect/
 
 tests:
 	# Running unittests
@@ -41,4 +41,9 @@ all:
 	@make html
 
 	# Running unittests
-	@make tests
+	@make tests
+
+
+
+develop:
+	pip install -e .

benchmarks/benchmark_cw.py

@@ -1,7 +1,8 @@
 """Benchmarks for CW complex computations"""
 import numpy as np
 import time
-from ect import ECT, EmbeddedCW, create_example_cw
+from ect import ECT
+from ect.utils.examples import create_example_cw
 import json
 from pathlib import Path
 
@@ -11,9 +12,9 @@ def benchmark_cw_ect(num_runs=5):
     results = {}
 
     configs = [
-        (8, 10),    # Small
-        (36, 36),   # Medium
-        (360, 360),  # Large
+        (8, 10),    
+        (36, 36),   
+        (360, 360),  
     ]
 
     for num_dir, num_thresh in configs:
@@ -26,7 +27,7 @@ def benchmark_cw_ect(num_runs=5):
             start_time = time.time()
 
             myect = ECT(num_dirs=num_dir, num_thresh=num_thresh)
-            myect.calculateECT(K)
+            myect.calculate(K)
 
             execution_time = time.time() - start_time
             times.append(execution_time)
@@ -45,7 +46,7 @@ def benchmark_cw_ect(num_runs=5):
     print("Running CW complex benchmarks...")
     results = benchmark_cw_ect()
 
-    # Save results
+
     output_dir = Path("benchmark_results")
     output_dir.mkdir(exist_ok=True)
 

benchmarks/benchmark_graph.py

@@ -21,12 +21,20 @@ def benchmark_graph_ect(num_runs=5):
     """Benchmark ECT computation on graphs"""
     results = {}
 
+    # Warmup run to trigger JIT compilation
+    warmup_shape = create_test_shape(100, 1)
+    G_warmup = EmbeddedGraph()
+    G_warmup.add_cycle(warmup_shape)
+    myect = ECT(num_dirs=360, num_thresh=360)
+    myect.calculate(G_warmup)  # Warmup run
+
     configs = [
         (100, 1),
         (1000, 1),
         (100, 3),
         (1000, 3),
         (10000, 3),
+        (100000, 5),
     ]
 
     for points, complexity in configs:
@@ -41,7 +49,7 @@ def benchmark_graph_ect(num_runs=5):
         for _ in range(num_runs):
             start_time = time.time()
             myect = ECT(num_dirs=360, num_thresh=360)
-            myect.calculateECT(G)
+            myect.calculate(G)
             times.append(time.time() - start_time)
 
         results[f'points_{points}_complexity_{complexity}'] = {
@@ -50,30 +58,9 @@ def benchmark_graph_ect(num_runs=5):
             'min_time': float(np.min(times)),
             'max_time': float(np.max(times))
         }
+    print(results)
 
     return results
 
 
-def benchmark_g_omega(num_runs=5):
-    """Benchmark g_omega computation"""
-    results = {}
-
-    sizes = [100, 500, 1000]
-    for size in sizes:
-        shape = create_test_shape(size)
-        G = EmbeddedGraph()
-        G.add_cycle(shape)
-
-        times = []
-        for _ in range(num_runs):
-            start_time = time.time()
-            for theta in np.linspace(0, 2*np.pi, 360):
-                G.g_omega(theta)
-            times.append(time.time() - start_time)
-
-        results[f'size_{size}'] = {
-            'mean_time': float(np.mean(times)),
-            'std_time': float(np.std(times))
-        }
 
-    return results

benchmarks/run_benchmarks.py

@@ -3,7 +3,7 @@
 import time
 from pathlib import Path
 import json
-from benchmark_graph import benchmark_graph_ect, benchmark_g_omega
+from benchmark_graph import benchmark_graph_ect
 from benchmark_cw import benchmark_cw_ect
 import platform
 
@@ -26,8 +26,7 @@ def run_all_benchmarks(num_runs=5):
     print("\nRunning CW complex benchmarks...")
     results['benchmarks']['cw_ect'] = benchmark_cw_ect(num_runs=num_runs)
 
-    print("\nRunning g_omega benchmarks...")
-    results['benchmarks']['g_omega'] = benchmark_g_omega(num_runs=num_runs)
+
 
     return results
 

doc_source/notebooks/tutorial_cw.py

@@ -0,0 +1,74 @@
+# %% [markdown]
+# # Tutorial: ECT for CW complexes
+#
+# This tutorial walks you through how to build a CW complex with the `EmbeddedCW` class, and then use the `ECT` class to compute the Euler characteristic transform
+
+# %%
+from ect import ECT, EmbeddedCW
+from ect.utils.examples import create_example_cw
+import numpy as np
+
+# %% [markdown]
+# The CW complex is the same as the `EmbeddedGraph` class with that additional ability to add faces. Faces are added by passing in a list of vertices. Note that we are generally assuming that these vertices follow around an empty region (as in, no other vertex is in the interior) in the graph bounded by the vertices, and further that all edges are already included in the graph. However the class does not yet check for this so you need to be careful!
+
+# %%
+K = EmbeddedCW()
+
+# Add vertices with coordinates
+K.add_node("A", [0, 0])
+K.add_node("B", [1, 0])
+K.add_node("C", [1, 1])
+K.add_node("D", [0, 1])
+
+# Add edges to form a square
+K.add_edges_from([("A", "B"), ("B", "C"), ("C", "D"), ("D", "A")])
+
+# Add the square face
+K.add_face(["A", "B", "C", "D"])
+
+K.center_coordinates()
+K.plot()
+
+# %% [markdown]
+# Just to have something a bit more interesting, let's make a more complicated example that's built into the class.
+
+# %%
+K = create_example_cw()
+K.plot(bounding_circle=True)
+
+
+# %% [markdown]
+# As with the `EmbeddedGraph` class, we can initialize the `ECT` class by deciding how many directions and how many thresholds to use.
+
+# %%
+ect = ECT(num_dirs=100, num_thresh=80)
+
+# %% [markdown]
+# Then we can compute the ECC for a single direction. In this case, the $x$-axis will be computed for the `num_thresh=80` stopping points in the interval $[-1.2r,1.2r]$ where $r$ is the minimum bounding radius for the input complex.
+
+# %%
+result = ect.calculate(K, theta=0)
+result.plot()
+
+# %% [markdown]
+# But of course it's easier to see this in a plot. This command calculates the ECC and immediately plots it.
+
+# %%
+# Set bounding radius to control threshold range
+radius = K.get_bounding_radius()
+ect.set_bounding_radius(1.2 * radius)
+
+# %% [markdown]
+# Similarly, we can compute the ECT and return the matrix. We make sure to internally set the bounding radius to use to control the $y$ axis of the plot.
+
+# %%
+result = ect.calculate(K)
+result.plot()
+
+# %% [markdown]
+# We can also look at the Smooth ECT:
+
+# %%
+# Calculate SECT and plot
+smooth = result.smooth()
+smooth.plot()