"""Boundary-clipped sunflower layout on the Fermat spiral of Vogel (1979).
Reference: https://doi.org/10.1016/0025-5564(79)90080-4
"""
from __future__ import annotations
from typing import TYPE_CHECKING
import numpy as np
from .._validation import _validate_finite
from .candidate import _CandidateLayout, _reach
from .fitters import SelectedCandidateFitter
if TYPE_CHECKING: # pragma: no cover
from ..boundaries.base import Boundary
__all__ = ["Sunflower"]
#: Golden angle in radians, ``pi (3 - sqrt 5)``.
_GOLDEN_ANGLE = float(np.pi * (3.0 - np.sqrt(5.0)))
[docs]
class Sunflower(_CandidateLayout):
r"""Fill a boundary with an exact-count phyllotactic pattern.
Candidates follow the Fermat spiral
.. math:: r_i=c\sqrt{i},\qquad \theta_i=i\gamma+\theta_0,
where ``c`` is radial spacing and ``gamma`` is the golden angle. The spiral
is centered on an interior point, which always holds a turbine, and clipped
by the boundary. The fitter finds the largest ``c`` leaving at least
``n_turbines`` candidates, then keeps an exact, well-distributed subset.
Parameters
----------
n_turbines:
Exact number of turbines to generate.
density:
The one density setting the count instead of ``n_turbines``, such as
``{"minimum_spacing": meters}``, as in :class:`LayoutAlgorithm`.
rotation:
Rotation of the spiral in degrees.
fitter:
Candidate-field fitting strategy. Defaults to one deterministic phase.
seed:
Random seed used only to resolve equivalent exact-count subsets.
"""
refine_candidates = False
def __init__(
self,
n_turbines: int | None = None,
*,
density: dict[str, float] | None = None,
rotation: float = 0.0,
fitter: SelectedCandidateFitter | None = None,
seed: int | np.random.Generator | None = None,
):
self.rotation = _validate_finite("rotation", rotation)
if fitter is None:
fitter = SelectedCandidateFitter(n_phases=1)
super().__init__(
n_turbines,
density=density,
fitter=fitter,
seed=seed,
)
@staticmethod
def _center_point(boundary) -> np.ndarray:
return np.asarray(boundary.interior_point, dtype=float)
[docs]
def estimate_spacing(self, boundary: Boundary) -> float:
"""Return a starting radial spacing ``c`` in meters for the fitter.
Assumes one ``area / n_turbines`` disc per turbine.
"""
return float(np.sqrt(boundary.area / (np.pi * self._required_n)))
def _candidate_origin(self, boundary):
return self._center_point(boundary)
def _required_candidates(self, boundary):
return self._center_point(boundary).reshape(1, 2)
def _candidate_bounds(self, boundary, scale, phases, max_candidates):
# Point i lies at radius scale * sqrt(i), and none past the boundary's
# reach from the center lies inside, so a field holds at most the
# points it builds, which the fitter's scan passes over where too few.
reach = _reach(boundary, self._center_point(boundary))
def bound(scales):
scales = np.asarray(scales, dtype=float)
# Capped below 2**63, so a tiny scale stays a large integer.
counts = np.minimum(np.ceil((reach / scales) ** 2) + 1, 2.0**62)
return np.repeat(counts.astype(np.int64)[:, None], len(phases), axis=1)
return bound
def _candidate_fields(
self, boundary, scale, phases, max_candidates=None, min_candidates=None
):
# The spiral has no phase, so every phase shares one field.
if scale <= 0:
return [np.empty((0, 2)) for _ in phases]
center = self._center_point(boundary)
maximum_radius = _reach(boundary, center)
maximum_index = int(np.ceil((maximum_radius / scale) ** 2))
if max_candidates is not None and maximum_index + 1 > max_candidates:
return [None for _ in phases]
indices = np.arange(maximum_index + 1, dtype=float)
angles = indices * _GOLDEN_ANGLE + np.deg2rad(self.rotation)
radii = scale * np.sqrt(indices)
points = center + np.column_stack(
[radii * np.cos(angles), radii * np.sin(angles)]
)
field = points[boundary._contains(points)]
return [field for _ in phases]
def _detailed_fit(self, boundary, positions, fit):
# The spiral has no phase, and its record names the pattern instead.
fit = dict(fit)
fit.pop("phase", None)
fit["method"] = "sunflower"
fit["pattern"] = "sunflower"
fit["rotation"] = self.rotation
fit["golden_angle_degrees"] = float(np.rad2deg(_GOLDEN_ANGLE))
return positions, fit
