moving the brownian curve implementation + associated refactorings

geosss/sphere.py

@@ -128,82 +128,3 @@ def rotate(theta, u=u, v=v, x=x):
         )
 
     return rotate
-
-
-def brownian_curve(n_points=100, dimension=6, step_size=0.05, seed=1234):
-    """
-    Generate smooth points on a unit d-sphere using Brownian motion.
-
-    Parameters:
-    n_points: number of points to generate
-    dimension: dimension of the sphere
-    step_size: step size for the Brownian motion, increasing it will make the points jump more
-    seed: random seed
-    """
-    rng = np.random.default_rng(seed)
-
-    # Initialize the first point on the dimension-1 -sphere
-    points = np.zeros((n_points, dimension))
-    points[0] = radial_projection(rng.standard_normal(dimension))
-
-    # Generate subsequent points via Brownian motion (small random steps)
-    for i in range(1, n_points):
-        step = rng.normal(size=dimension) * step_size
-        new_point = points[i - 1] + step
-
-        # Project back to the unit sphere
-        points[i] = radial_projection(new_point)
-
-    return points
-
-
-def constrained_brownian_curve(n_points=100, dimension=6, step_size=0.05, seed=1234):
-    """
-    Generate smooth points on a unit d-sphere using constrained Brownian motion
-    to avoid loops and overlaps.
-
-    Parameters:
-    n_points: number of points to generate
-    dimension: dimension of the sphere (must be at least 2)
-    step_size: step size for the motion
-    seed: random seed
-    """
-    rng = np.random.default_rng(seed)
-
-    # Initialize the first point on the (dimension - 1)-sphere
-    points = np.zeros((n_points, dimension))
-    points[0] = radial_projection(rng.standard_normal(dimension))
-
-    # Initialize the direction of motion (tangent vector at the first point)
-    v = rng.standard_normal(dimension)
-    v -= np.dot(v, points[0]) * points[0]  # Make v orthogonal to points[0]
-    v /= np.linalg.norm(v)  # Normalize the tangent vector
-
-    for i in range(1, n_points):
-        # Generate a small random perturbation orthogonal to both v and points[i - 1]
-        random_step = rng.standard_normal(dimension)
-
-        # Make the random perturbation tangent to the sphere at points[i - 1]
-        random_step -= np.dot(random_step, points[i - 1]) * points[i - 1]
-
-        # Make the random perturbation orthogonal to the current direction v
-        random_step -= np.dot(random_step, v) * v
-
-        # Normalize the random_step so that it lies on the tangent plane
-        random_step /= np.linalg.norm(random_step)
-
-        # Scale the random_step with the given step size
-        random_step *= step_size
-
-        # Update the tangent vector v
-        v += random_step
-
-        # Ensure the updated tangent vector v remains tangent to the sphere
-        v -= np.dot(v, points[i - 1]) * points[i - 1]
-        v /= np.linalg.norm(v)
-
-        # Move the point along v by a small angle (step_size)
-        angle = step_size
-        points[i] = points[i - 1] * np.cos(angle) + v * np.sin(angle)
-
-    return points

geosss/spherical_curve.py

@@ -100,3 +100,82 @@ def find_nearest(self, point):
             distances.append(d)
             points.append(y)
         return points[np.argmin(distances)]
+
+
+def brownian_curve(n_points=100, dimension=6, step_size=0.05, seed=1234):
+    """
+    Generate smooth points on a unit d-sphere using Brownian motion.
+
+    Parameters:
+    n_points: number of points to generate
+    dimension: dimension of the sphere
+    step_size: step size for the Brownian motion, increasing it will make the points jump more
+    seed: random seed
+    """
+    rng = np.random.default_rng(seed)
+
+    # Initialize the first point on the dimension-1 -sphere
+    points = np.zeros((n_points, dimension))
+    points[0] = sphere.radial_projection(rng.standard_normal(dimension))
+
+    # Generate subsequent points via Brownian motion (small random steps)
+    for i in range(1, n_points):
+        step = rng.normal(size=dimension) * step_size
+        new_point = points[i - 1] + step
+
+        # Project back to the unit sphere
+        points[i] = sphere.radial_projection(new_point)
+
+    return points
+
+
+def constrained_brownian_curve(n_points=100, dimension=6, step_size=0.05, seed=1234):
+    """
+    Generate smooth points on a unit d-sphere using constrained Brownian motion
+    to avoid loops and overlaps.
+
+    Parameters:
+    n_points: number of points to generate
+    dimension: dimension of the sphere (must be at least 2)
+    step_size: step size for the motion
+    seed: random seed
+    """
+    rng = np.random.default_rng(seed)
+
+    # Initialize the first point on the (dimension - 1)-sphere
+    points = np.zeros((n_points, dimension))
+    points[0] = sphere.radial_projection(rng.standard_normal(dimension))
+
+    # Initialize the direction of motion (tangent vector at the first point)
+    v = rng.standard_normal(dimension)
+    v -= np.dot(v, points[0]) * points[0]  # Make v orthogonal to points[0]
+    v /= np.linalg.norm(v)  # Normalize the tangent vector
+
+    for i in range(1, n_points):
+        # Generate a small random perturbation orthogonal to both v and points[i - 1]
+        random_step = rng.standard_normal(dimension)
+
+        # Make the random perturbation tangent to the sphere at points[i - 1]
+        random_step -= np.dot(random_step, points[i - 1]) * points[i - 1]
+
+        # Make the random perturbation orthogonal to the current direction v
+        random_step -= np.dot(random_step, v) * v
+
+        # Normalize the random_step so that it lies on the tangent plane
+        random_step /= np.linalg.norm(random_step)
+
+        # Scale the random_step with the given step size
+        random_step *= step_size
+
+        # Update the tangent vector v
+        v += random_step
+
+        # Ensure the updated tangent vector v remains tangent to the sphere
+        v -= np.dot(v, points[i - 1]) * points[i - 1]
+        v /= np.linalg.norm(v)
+
+        # Move the point along v by a small angle (step_size)
+        angle = step_size
+        points[i] = points[i - 1] * np.cos(angle) + v * np.sin(angle)
+
+    return points

scripts/curve_Nd.py

@@ -10,7 +10,7 @@
 
 import geosss as gs
 from geosss.distributions import CurvedVonMisesFisher, Distribution
-from geosss.spherical_curve import SlerpCurve, SphericalCurve
+from geosss.spherical_curve import SlerpCurve, SphericalCurve, brownian_curve
 
 mpl.rcParams["mathtext.fontset"] = "cm"  # Use Computer Modern font
 
@@ -307,7 +307,7 @@ def argparser():
     n_runs = args.n_runs  # sampler runs (default: 1), for ess computations `n_runs=10`
 
     # optional controls
-    brownian_curve = True  # brownian curve or curve with fixed knots
+    is_brownian_curve = True  # brownian curve or curve with fixed knots
     reprod_switch = True  # seeds samplers for reproducibility
     show_plots = True  # show the plots
     savefig = True  # save the plots
@@ -319,7 +319,7 @@ def argparser():
     setup_logging(savedir, kappa)
 
     # Define curve on the sphere
-    if not brownian_curve:
+    if not is_brownian_curve:
         knots = np.array(
             [
                 [-0.25882694, 0.95006168, 0.17433133],
@@ -339,18 +339,16 @@ def argparser():
         if n_dim > knots.shape[1]:
             knots = np.pad(knots, ((0, 0), (n_dim - knots.shape[1], 0)))
 
-        curve = SlerpCurve(knots)
-
     else:
         # generates a smooth curve on the sphere with brownian motion
-        knots = gs.sphere.brownian_curve(
+        knots = brownian_curve(
             n_points=10,
             dimension=n_dim,
             step_size=0.5,  # larger step size will result in more spread out points
             seed=4562,
         )
-        curve = SlerpCurve(knots)
-        # curve = SlerpCurve.random_curve(n_knots=100, seed=4562, dimension=n_dim)
+
+    curve = SlerpCurve(knots)
 
     # Initialize based on dimensionality
     initial = (

scripts/curve_Nd_parallelized_nruns.py

@@ -8,7 +8,7 @@
 
 import geosss as gs
 from geosss.distributions import CurvedVonMisesFisher, Distribution
-from geosss.spherical_curve import SlerpCurve
+from geosss.spherical_curve import SlerpCurve, brownian_curve
 
 
 def setup_logging(
@@ -298,12 +298,12 @@ def main():
     assert n_dim >= 3, msg
 
     # optional controls
-    brownian_curve = True  # brownian curve or curve with fixed knots
+    is_brownian_curve = True  # brownian curve or curve with fixed knots
     reprod_switch = True  # seeds samplers for reproducibility
     rerun_if_samples_exists = False  # rerun even if samples file exists
 
     # creating a target as a curve on the sphere
-    if not brownian_curve:
+    if not is_brownian_curve:
         knots = np.array(
             [
                 [-0.25882694, 0.95006168, 0.17433133],
@@ -325,13 +325,14 @@ def main():
 
     else:
         # generates a smooth curve on the sphere with brownian motion
-        knots = gs.sphere.brownian_curve(
+        knots = brownian_curve(
             n_points=10,
             dimension=n_dim,
             step_size=0.5,  # larger step size will result in more spread out points
             seed=4562,
         )
 
+    logging.info(f"Target curve: {knots}")
     # defining the curve on the sphere
     curve = SlerpCurve(knots)