Fractals

Fractals: escape-time sets, chaos-game attractors, and words and patches of named curves and sets.

CantorElement(label, tiles=None, *, ratio=1 / 3, prototiles=(), filled=False, max_size=10000000, bounds=None, visual=None)

Bases: PatchElement

Patch of intervals, each replaced by its first and its last part of the length ratio.

Parameters:

Name	Type	Description	Default
label	str	The label of the patch.	required
tiles	Sequence[tuple[str, AffineSpec]] | None	Pairs of a prototile name and its placement at generation 0, or None for the first prototile in place.	None
ratio	Annotated[float, ParameterSpec(domain=open(0, oo))]	The length of each kept part relative to its interval.	1 / 3
prototiles	Sequence[PrototileSpec]	The prototiles of the patch.	()
filled	bool	Whether the geometry visual draws polygons as a filled image instead of outlines.	False
max_size	int	The largest number of tiles a generation may produce.	10000000
bounds	Sequence[float] | None	The box (xmin, ymin, xmax, ymax) the geometry visual shows, or None for the box of the initial tiles.	None
visual	VisualSpec | None	How the patch is drawn, or None for default_visual_.	None

Examples:

>>> from skcomplex.systems.mathematics.fractals import CantorElement
>>> CantorElement('cantor', ratio=0.25).ratio
0.25

Source code in src/skcomplex/systems/mathematics/fractals/_substitution.py

def __init__(
    self,
    label: str,
    tiles: Sequence[tuple[str, AffineSpec]] | None = None,
    *,
    ratio: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1 / 3,
    prototiles: Sequence[PrototileSpec] = (),
    filled: bool = False,
    max_size: int = 10_000_000,
    bounds: Sequence[float] | None = None,
    visual: VisualSpec | None = None,
) -> None:
    super().__init__(
        label,
        tiles,
        prototiles=prototiles,
        filled=filled,
        max_size=max_size,
        bounds=bounds,
        visual=visual,
    )
    self.ratio = ratio

ChaosPointElement(label, visual=VisualSpec('points', {'position': 'p', 'radius': 0.004, 'states': ()}), *, p=0.0)

Bases: Element

Point of the plane whose position a random affine map moves at every step.

Parameters:

Name	Type	Description	Default
label	str	The label of the point or the group.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	VisualSpec('points', {'position': 'p', 'radius': 0.004, 'states': ()})
p	Annotated[object, StateSpec(vector=True)]	The initial position.	0.0

Examples:

>>> from skcomplex.systems.mathematics.fractals import ChaosPointElement
>>> ChaosPointElement('point', p=(0.5, 0.5)).p
(0.5, 0.5)

Source code in src/skcomplex/systems/mathematics/fractals/_chaos_game.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('points', {'position': 'p', 'radius': 0.004, 'states': ()}),
    *,
    p: Annotated[object, StateSpec(vector=True)] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.p = p

InverseJuliaPointElement(label, visual=None, *, z=0.5, omega=0.0, p=0.5)

Bases: Element

Point of the inverse iteration of z -> z² + omega, which moves to a square root of z - omega drawn at random.

The point moves to +√(z - omega) with probability p and to -√(z - omega) otherwise, with the principal square root.

Parameters:

Name	Type	Description	Default
label	str	The label of the element or the group.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
z	Annotated[object, StateSpec(domain=CC)]	The initial point.	0.5
omega	Annotated[object, ParameterSpec(domain=CC)]	The constant of the quadratic map.	0.0
p	Annotated[object, ParameterSpec(domain=closed(0, 1))]	The probability of the positive square root.	0.5

Examples:

>>> from skcomplex.systems.mathematics.fractals import InverseJuliaPointElement
>>> InverseJuliaPointElement('point', p=0.3).p
0.3

Source code in src/skcomplex/systems/mathematics/fractals/_quadratic.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    z: Annotated[object, StateSpec(domain=CC)] = 0.5,
    omega: Annotated[object, ParameterSpec(domain=CC)] = 0.0,
    p: Annotated[object, ParameterSpec(domain=RealSet.closed(0, 1))] = 0.5,
) -> None:
    super().__init__(label, visual=visual)
    self.z = z
    self.omega = omega
    self.p = p

KochElement(label, axiom='', *, alphabet=None, ignore='', geometry=None, max_size=10000000, bounds=None, visual=None)

Bases: WordElement

Word rewritten by the production F → F+F--F+F of the Koch curve.

Parameters:

Name	Type	Description	Default
label	str	The label of the word.	required
axiom	Annotated[Letters, StateSpec(domain=None)]	The word at generation 0.	''
geometry	GeometrySpec | None	How the word is drawn as line segments, or None for no geometry.	None

Examples:

>>> from skcomplex.systems.mathematics.fractals import KochElement
>>> KochElement('koch', axiom='F')._get_laws(None)[0].successor
'F+F--F+F'

Source code in src/skcomplex/core/components/elements/_word.py

def __init__(
    self,
    label: str,
    axiom: Annotated[Letters, StateSpec(domain=None)] = '',
    *,
    alphabet: str | Mapping[str, int] | None = None,
    ignore: str | Sequence[str] = '',
    geometry: GeometrySpec | None = None,
    max_size: int = 10_000_000,
    bounds: Sequence[float] | None = None,
    visual: VisualSpec | None = None,
) -> None:
    super().__init__(label, visual=visual)
    self.axiom = axiom
    self.alphabet = alphabet
    self.ignore = ignore
    self.geometry = geometry
    self.max_size = max_size
    self.bounds = bounds

QuadraticMapElement(label, visual=None, *, z=0.0, n=0.0, c=0.0, radius=2.0)

Bases: Element

Complex quadratic map z' = z^2 + c while |z| stays within the escape radius, with a count of its steps.

Parameters:

Name	Type	Description	Default
label	str	The label of the element.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
z	Annotated[object, StateSpec(domain=CC)]	The initial point.	0.0
n	Annotated[object, StateSpec(domain=NN)]	The initial number of steps within the escape radius.	0.0
c	Annotated[object, ParameterSpec(domain=CC)]	The constant of the map.	0.0
radius	Annotated[object, ParameterSpec(domain=open(0, oo))]	The escape radius.	2.0

Examples:

>>> from skcomplex.systems.mathematics.fractals import QuadraticMapElement
>>> QuadraticMapElement('point', c=0.25).c
0.25

Source code in src/skcomplex/systems/mathematics/fractals/_quadratic.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    z: Annotated[object, StateSpec(domain=CC)] = 0.0,
    n: Annotated[object, StateSpec(domain=NN)] = 0.0,
    c: Annotated[object, ParameterSpec(domain=CC)] = 0.0,
    radius: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 2.0,
) -> None:
    super().__init__(label, visual=visual)
    self.z = z
    self.n = n
    self.c = c
    self.radius = radius

RandomAffineMapInteraction(label, *, target, state='p', maps=SIERPINSKI_MAPS, weights=(1.0, 1.0, 1.0))

Bases: MapInteraction

Map of a point of the plane by one affine map of an iterated function system, drawn anew at every step.

Map k takes (x, y) to (a x + b y + e, c x + d y + f) for its coefficients (a, b, c, d, e, f), and is drawn with a probability proportional to its weight, independently for every member of a bound group.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction, unique among the parts of its system.	required
target	Annotated[str, RoleSpec()]	The path of the element or group whose position the map sets.	required
state	str	The name of the vector state of the target that holds the position.	'p'
maps	Sequence[AffineMap]	The coefficients (a, b, c, d, e, f) of every affine map.	SIERPINSKI_MAPS
weights	Sequence[float]	The weight of every map.	(1.0, 1.0, 1.0)

Raises:

Type	Description
ValidationError	The number of weights differs from the number of maps, or a map does not have six coefficients.

Examples:

>>> from skcomplex.systems.mathematics.fractals import RandomAffineMapInteraction
>>> game = RandomAffineMapInteraction('game', target='points', maps=[(0.5, 0, 0, 0.5, 0, 0)], weights=[1.0])
>>> game.maps
[(0.5, 0, 0, 0.5, 0, 0)]

Source code in src/skcomplex/systems/mathematics/fractals/_chaos_game.py

def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    state: str = 'p',
    maps: Sequence[AffineMap] = SIERPINSKI_MAPS,
    weights: Sequence[float] = (1.0, 1.0, 1.0),
) -> None:
    super().__init__(label, target=target, state=state, function=self._get_image)
    self.maps = maps
    self.weights = weights

SierpinskiElement(label, tiles=None, *, shape='triangle', prototiles=(), filled=False, max_size=10000000, bounds=None, visual=None)

Bases: PatchElement

Patch of triangles or squares, each replaced by the scaled copies of the Sierpiński triangle or carpet.

Parameters:

Name	Type	Description	Default
label	str	The label of the patch.	required
tiles	Sequence[tuple[str, AffineSpec]] | None	Pairs of a prototile name and its placement at generation 0, or None for the first prototile in place.	None
shape	str	'triangle' or 'carpet'.	'triangle'
prototiles	Sequence[PrototileSpec]	The prototiles of the patch.	()
filled	bool	Whether the geometry visual draws polygons as a filled image instead of outlines.	False
max_size	int	The largest number of tiles a generation may produce.	10000000
bounds	Sequence[float] | None	The box (xmin, ymin, xmax, ymax) the geometry visual shows, or None for the box of the initial tiles.	None
visual	VisualSpec | None	How the patch is drawn, or None for default_visual_.	None

Examples:

>>> from skcomplex.systems.mathematics.fractals import SierpinskiElement
>>> len(SierpinskiElement('sierpinski', shape='carpet')._get_laws(None)[0].pieces)
8

Source code in src/skcomplex/systems/mathematics/fractals/_substitution.py

def __init__(
    self,
    label: str,
    tiles: Sequence[tuple[str, AffineSpec]] | None = None,
    *,
    shape: str = 'triangle',
    prototiles: Sequence[PrototileSpec] = (),
    filled: bool = False,
    max_size: int = 10_000_000,
    bounds: Sequence[float] | None = None,
    visual: VisualSpec | None = None,
) -> None:
    super().__init__(
        label,
        tiles,
        prototiles=prototiles,
        filled=filled,
        max_size=max_size,
        bounds=bounds,
        visual=visual,
    )
    self.shape = shape

cantor_set(ratio=1 / 3)

Return the Cantor set as the substitution of every interval by its two outer parts.

Every generation keeps the first and the last part of each interval, each ratio of its length, and removes the middle. With the ratio of one third, generation n has 2^n intervals of total length (2/3)^n, and the limit set has the similarity dimension log 2 / log 3, about 0.631.

Parameters:

Name	Type	Description	Default
ratio	float	The length of each kept part relative to its interval, below one half.	1 / 3

Returns:

Name	Type	Description
system	System	The system, with the patch cantor and its number of intervals cantor.size.

Raises:

Type	Description
ValidationError	The ratio is not between 0 and 1/2.

Examples:

>>> from skcomplex.systems import create_simulation
>>> create_simulation('cantor_set').prepare().run(until=3).trajectory_.values('cantor.size').tolist()
[1.0, 2.0, 4.0, 8.0]

Source code in src/skcomplex/systems/mathematics/fractals/_substitution.py

@register_system(
    field='mathematics/fractals',
    tags=('fractal', 'substitution', 'rewriting system'),
    references=('Cantor 1883', 'Thurner, Hanel, and Klimek 2018', 'Trzaska 2023'),
    defaults={'compiler': RewritingCompiler(), 'until': 8, 'dt': 1},
)
def cantor_set(ratio: float = 1 / 3) -> System:
    """Return the Cantor set as the substitution of every interval by its two outer parts.

    Every generation keeps the first and the last part of each interval, each `ratio` of its length, and removes the
    middle. With the ratio of one third, generation `n` has `2^n` intervals of total length `(2/3)^n`, and the limit set
    has the similarity dimension `log 2 / log 3`, about 0.631.

    Args:
        ratio:
            The length of each kept part relative to its interval, below one half.

    Returns:
        system:
            The system, with the patch `cantor` and its number of intervals `cantor.size`.

    Raises:
        ValidationError:
            The ratio is not between 0 and 1/2.

    Examples:
        >>> from skcomplex.systems import create_simulation
        >>> create_simulation('cantor_set').prepare().run(until=3).trajectory_.values('cantor.size').tolist()
        [1.0, 2.0, 4.0, 8.0]
    """
    if not 0 < ratio < LARGEST_RATIO:
        error_msg = f'ratio = {ratio!r} must lie between 0 and 1/2.'
        raise ValidationError(error_msg)
    interval = PrototileSpec('interval', [(0.0, 0.0), (1.0, 0.0)])
    return System(
        elements=[CantorElement('cantor', prototiles=[interval], bounds=(-0.03, -0.1, 1.03, 0.1), ratio=ratio)],
    )

chaos_game(preset='sierpinski', *, maps=None, weights=None, n_points=20000, random_state=0)

Return points that draw the attractor of an iterated function system by the chaos game.

At every step each point moves by one of the affine maps of the system, drawn with probabilities proportional to the weights. The maps contract, so after a few steps the cloud of points covers the attractor, the unique compact set that the union of the maps leaves unchanged. The Sierpiński triangle of three maps of ratio 1/2 has the box-counting dimension log 3 / log 2, and the Barnsley fern is the attractor of four maps.

Parameters:

Name	Type	Description	Default
preset	str	'sierpinski' or 'barnsley', used when maps is None.	'sierpinski'
maps	Sequence[AffineMap] | None	The coefficients (a, b, c, d, e, f) of every map (x, y) -> (a x + b y + e, c x + d y + f), or None for the preset.	None
weights	Sequence[float] | None	The weight of every map, or None for the weights of the preset, or equal weights for given maps.	None
n_points	int	The number of points, started uniformly in the unit square.	20000
random_state	int | Generator | None	The seed or generator of the initial points.	0

Returns:

Name	Type	Description
system	System	The system, with the positions points.p.

Raises:

Type	Description
ValidationError	The preset is unknown.

Examples:

>>> from skcomplex.systems.mathematics.fractals import chaos_game
>>> len(chaos_game('barnsley', n_points=10).interactions[0].maps)
4

Source code in src/skcomplex/systems/mathematics/fractals/_chaos_game.py

@register_system(
    field='mathematics/fractals',
    tags=('fractal', 'iterated function system', 'chaos game', 'discrete map'),
    references=('Barnsley 1988', 'Dilão 2023'),
    defaults={'until': 30, 'dt': 1},
)
def chaos_game(
    preset: str = 'sierpinski',
    *,
    maps: Sequence[AffineMap] | None = None,
    weights: Sequence[float] | None = None,
    n_points: int = 20000,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return points that draw the attractor of an iterated function system by the chaos game.

    At every step each point moves by one of the affine maps of the system, drawn with probabilities proportional to
    the weights. The maps contract, so after a few steps the cloud of points covers the attractor, the unique compact
    set that the union of the maps leaves unchanged. The Sierpiński triangle of three maps of ratio 1/2 has the
    box-counting dimension `log 3 / log 2`, and the Barnsley fern is the attractor of four maps.

    Args:
        preset:
            `'sierpinski'` or `'barnsley'`, used when `maps` is `None`.

        maps:
            The coefficients `(a, b, c, d, e, f)` of every map `(x, y) -> (a x + b y + e, c x + d y + f)`, or `None` for
            the preset.

        weights:
            The weight of every map, or `None` for the weights of the preset, or equal weights for given maps.

        n_points:
            The number of points, started uniformly in the unit square.

        random_state:
            The seed or generator of the initial points.

    Returns:
        system:
            The system, with the positions `points.p`.

    Raises:
        ValidationError:
            The preset is unknown.

    Examples:
        >>> from skcomplex.systems.mathematics.fractals import chaos_game
        >>> len(chaos_game('barnsley', n_points=10).interactions[0].maps)
        4
    """
    if maps is None and preset not in PRESETS:
        error_msg = f'Unknown chaos game preset {preset!r}. Use one of {", ".join(PRESETS)}, or give the maps.'
        raise ValidationError(error_msg)
    chosen_maps, preset_weights = PRESETS[preset] if maps is None else (maps, (1.0,) * len(maps))
    return _system_of(chosen_maps, weights if weights is not None else preset_weights, n_points, random_state)

inverse_julia(omega=-0.123 + 0.745j, *, p=0.5, n_points=5000, random_state=0)

Return the random function system z -> ±√(z - omega), whose points settle on the Julia set of z² + omega.

At every step each point moves to the positive square root of z - omega with probability p and to the negative one otherwise (Dilão 2023, Exercise 5.8). Both branches invert the quadratic map z -> z² + omega, which repels its Julia set, so the inverse branches attract the points onto the Julia set, drawn by the chaos game with weights p and 1 - p. For omega = 0 the Julia set is the unit circle, and the default is Douady's rabbit.

Parameters:

Name	Type	Description	Default
omega	complex	The constant of the quadratic map.	-0.123 + 0.745j
p	float	The probability of the positive square root.	0.5
n_points	int	The number of points, started uniformly in the square of side 2 about the origin.	5000
random_state	int | Generator | None	The seed or generator of the initial points.	0

Returns:

Name	Type	Description
system	System	The system, with the points points.z.

Examples:

>>> from skcomplex.systems.mathematics.fractals import inverse_julia
>>> inverse_julia(0.0, n_points=10).elements[0].prepare().size_
10

Source code in src/skcomplex/systems/mathematics/fractals/_quadratic.py

@register_system(
    field='mathematics/fractals',
    tags=('fractal', 'complex', 'iterated function system', 'random map', 'julia set'),
    references=('Barnsley 1988', 'Dilão 2023'),
    defaults={'until': 30, 'dt': 1},
)
def inverse_julia(
    omega: complex = -0.123 + 0.745j,
    *,
    p: float = 0.5,
    n_points: int = 5000,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return the random function system `z -> ±√(z - omega)`, whose points settle on the Julia set of `z² + omega`.

    At every step each point moves to the positive square root of `z - omega` with probability `p` and to the
    negative one otherwise (Dilão 2023, Exercise 5.8). Both branches invert the quadratic map `z -> z² + omega`, which
    repels its Julia set, so the inverse branches attract the points onto the Julia set, drawn by the chaos game with
    weights `p` and `1 - p`. For `omega = 0` the Julia set is the unit circle, and the default is Douady's rabbit.

    Args:
        omega:
            The constant of the quadratic map.

        p:
            The probability of the positive square root.

        n_points:
            The number of points, started uniformly in the square of side 2 about the origin.

        random_state:
            The seed or generator of the initial points.

    Returns:
        system:
            The system, with the points `points.z`.

    Examples:
        >>> from skcomplex.systems.mathematics.fractals import inverse_julia
        >>> inverse_julia(0.0, n_points=10).elements[0].prepare().size_
        10
    """
    rng = np.random.default_rng(random_state)
    start = rng.uniform(-1.0, 1.0, n_points) + 1j * rng.uniform(-1.0, 1.0, n_points)
    values = {'z': start, 'omega': complex(omega), 'p': p}
    return System(elements=[GroupElement('points', element=InverseJuliaPointElement, size=n_points, values=values)])

julia(c=-0.8 + 0.156j, shape=(240, 320), *, window=JULIA_WINDOW, radius=2.0)

Return the filled Julia set of z -> z^2 + c, the starting points whose orbit stays bounded (Dilão 2023 §5.5).

Every site holds one starting point z of the window and iterates the map with the same constant c, and n counts the steps the orbit stays within the escape radius. For c in the Mandelbrot set the filled Julia set is connected, and outside it the set breaks into dust.

Parameters:

Name	Type	Description	Default
c	complex	The constant of the map.	-0.8 + 0.156j
shape	Sequence[int]	The number of sites along the imaginary and the real axis.	(240, 320)
window	Window	The smallest and largest real part, then the smallest and largest imaginary part, of the starting points.	JULIA_WINDOW
radius	float	The escape radius, at least 2 and at least |c|, so escaped orbits never return.	2.0

Returns:

Name	Type	Description
system	System	The system, with the field plane, its orbits plane.z, and its escape times plane.n.

Examples:

>>> from skcomplex.systems.mathematics.fractals import julia
>>> julia(0.0, (3, 3), window=(-1, 1, -1, 1)).elements[0].values['z'][0].tolist()
[(-1+1j), 1j, (1+1j)]

Source code in src/skcomplex/systems/mathematics/fractals/_quadratic.py

@register_system(
    field='mathematics/fractals',
    tags=('fractal', 'complex', 'map', 'lattice', 'escape time'),
    references=('Dilão 2023',),
    defaults={'until': 100, 'dt': 1},
)
def julia(
    c: complex = -0.8 + 0.156j,
    shape: Sequence[int] = (240, 320),
    *,
    window: Window = JULIA_WINDOW,
    radius: float = 2.0,
) -> System:
    """Return the filled Julia set of `z -> z^2 + c`, the starting points whose orbit stays bounded (Dilão 2023 §5.5).

    Every site holds one starting point `z` of the window and iterates the map with the same constant `c`, and `n`
    counts the steps the orbit stays within the escape radius. For `c` in the Mandelbrot set the filled Julia set is
    connected, and outside it the set breaks into dust.

    Args:
        c:
            The constant of the map.

        shape:
            The number of sites along the imaginary and the real axis.

        window:
            The smallest and largest real part, then the smallest and largest imaginary part, of the starting points.

        radius:
            The escape radius, at least 2 and at least `|c|`, so escaped orbits never return.

    Returns:
        system:
            The system, with the field `plane`, its orbits `plane.z`, and its escape times `plane.n`.

    Examples:
        >>> from skcomplex.systems.mathematics.fractals import julia
        >>> julia(0.0, (3, 3), window=(-1, 1, -1, 1)).elements[0].values['z'][0].tolist()
        [(-1+1j), 1j, (1+1j)]
    """
    points = _plane(shape, window)
    return _escape_field(shape, points, complex(c), radius)

koch_curve(angle=60.0)

Return the Koch curve as the L-system F → F+F--F+F drawn by a turtle.

Every generation replaces each segment by four segments of a third of its length, with a bump of the turning angle in the middle. With the angle of 60 degrees the curve keeps its endpoints at (0, 0) and (1, 0), its length is (4/3)^n at generation n, and its similarity dimension is log 4 / log 3, about 1.262.

Parameters:

Name	Type	Description	Default
angle	float	The turning angle of the turtle in degrees.	60.0

Returns:

Name	Type	Description
system	System	The system, with the word koch, its symbol counts koch.F, and its size koch.size.

Examples:

>>> from skcomplex.systems import create_simulation
>>> create_simulation('koch_curve').prepare().run(until=3).trajectory_.values('koch.F').tolist()
[1.0, 4.0, 16.0, 64.0]

Source code in src/skcomplex/systems/mathematics/fractals/_substitution.py

@register_system(
    field='mathematics/fractals',
    tags=('fractal', 'L-system', 'rewriting system', 'turtle'),
    references=('Prusinkiewicz and Lindenmayer 1990', 'Thurner, Hanel, and Klimek 2018', 'Trzaska 2023'),
    defaults={'compiler': RewritingCompiler(), 'until': 5, 'dt': 1},
)
def koch_curve(angle: float = 60.0) -> System:
    """Return the Koch curve as the L-system `F → F+F--F+F` drawn by a turtle.

    Every generation replaces each segment by four segments of a third of its length, with a bump of the turning angle
    in the middle. With the angle of 60 degrees the curve keeps its endpoints at `(0, 0)` and `(1, 0)`, its length is
    `(4/3)^n` at generation `n`, and its similarity dimension is `log 4 / log 3`, about 1.262.

    Args:
        angle:
            The turning angle of the turtle in degrees.

    Returns:
        system:
            The system, with the word `koch`, its symbol counts `koch.F`, and its size `koch.size`.

    Examples:
        >>> from skcomplex.systems import create_simulation
        >>> create_simulation('koch_curve').prepare().run(until=3).trajectory_.values('koch.F').tolist()
        [1.0, 4.0, 16.0, 64.0]
    """
    word = KochElement(
        'koch',
        axiom='F',
        alphabet='F+-',
        geometry=GeometrySpec(angle=angle, step=1.0, contraction=1 / 3),
        bounds=(-0.03, -0.05, 1.03, 0.35),
    )
    return System(elements=[word])

mandelbrot(shape=(240, 300), *, window=MANDELBROT_WINDOW, radius=2.0)

Return the Mandelbrot set, the points c whose orbit of z -> z^2 + c from 0 stays bounded (Dilão 2023 §5.5).

Every site holds one point c of the window and iterates the map, and n counts the steps the orbit stays within the escape radius. The points of the set never escape, so their count grows with every step, while the count of an escaped point stops at its escape time and colors the bands around the set.

Parameters:

Name	Type	Description	Default
shape	Sequence[int]	The number of sites along the imaginary and the real axis.	(240, 300)
window	Window	The smallest and largest real part, then the smallest and largest imaginary part, of the points.	MANDELBROT_WINDOW
radius	float	The escape radius, 2 for the quadratic map.	2.0

Returns:

Name	Type	Description
system	System	The system, with the field plane, its orbits plane.z, and its escape times plane.n.

Examples:

>>> from skcomplex.systems.mathematics.fractals import mandelbrot
>>> mandelbrot((3, 3), window=(-1, 1, -1, 1)).elements[0].values['c'][1].tolist()
[(-1+0j), 0j, (1+0j)]

Source code in src/skcomplex/systems/mathematics/fractals/_quadratic.py

@register_system(
    field='mathematics/fractals',
    tags=('fractal', 'complex', 'map', 'lattice', 'escape time'),
    references=('Dilão 2023',),
    defaults={'until': 100, 'dt': 1},
)
def mandelbrot(
    shape: Sequence[int] = (240, 300),
    *,
    window: Window = MANDELBROT_WINDOW,
    radius: float = 2.0,
) -> System:
    """Return the Mandelbrot set, the points `c` whose orbit of `z -> z^2 + c` from 0 stays bounded (Dilão 2023 §5.5).

    Every site holds one point `c` of the window and iterates the map, and `n` counts the steps the orbit stays within
    the escape radius. The points of the set never escape, so their count grows with every step, while the count of an
    escaped point stops at its escape time and colors the bands around the set.

    Args:
        shape:
            The number of sites along the imaginary and the real axis.

        window:
            The smallest and largest real part, then the smallest and largest imaginary part, of the points.

        radius:
            The escape radius, 2 for the quadratic map.

    Returns:
        system:
            The system, with the field `plane`, its orbits `plane.z`, and its escape times `plane.n`.

    Examples:
        >>> from skcomplex.systems.mathematics.fractals import mandelbrot
        >>> mandelbrot((3, 3), window=(-1, 1, -1, 1)).elements[0].values['c'][1].tolist()
        [(-1+0j), 0j, (1+0j)]
    """
    points = _plane(shape, window)
    return _escape_field(shape, np.zeros_like(points), points, radius)

sierpinski_substitution(shape='triangle', *, filled=True)

Return the Sierpiński triangle or carpet as the substitution of a polygon by scaled copies of itself.

The triangle is replaced by its three corner copies of half its side, so generation n has 3^n triangles of total area (3/4)^n times the first, and the gasket has the similarity dimension log 3 / log 2, about 1.585. The square is replaced by the eight copies of a third of its side around its middle, so generation n has 8^n squares of area (8/9)^n, and the carpet has the dimension log 8 / log 3, about 1.893.

Parameters:

Name	Type	Description	Default
shape	str	'triangle' or 'carpet'.	'triangle'
filled	bool	Whether the geometry visual draws the polygons as a filled image instead of outlines.	True

Returns:

Name	Type	Description
system	System	The system, with the patch sierpinski and its number of polygons sierpinski.size.

Raises:

Type	Description
ValidationError	The shape is unknown.

Examples:

>>> from skcomplex.systems import create_simulation
>>> simulation = create_simulation('sierpinski_substitution').prepare().run(until=2)
>>> simulation.trajectory_.values('sierpinski.size').tolist()
[1.0, 3.0, 9.0]

Source code in src/skcomplex/systems/mathematics/fractals/_substitution.py

@register_system(
    field='mathematics/fractals',
    tags=('fractal', 'substitution', 'rewriting system'),
    references=('Sierpiński 1915', 'Sierpiński 1916', 'Thurner, Hanel, and Klimek 2018', 'Trzaska 2023'),
    defaults={'compiler': RewritingCompiler(), 'until': 6, 'dt': 1},
)
def sierpinski_substitution(shape: str = 'triangle', *, filled: bool = True) -> System:
    """Return the Sierpiński triangle or carpet as the substitution of a polygon by scaled copies of itself.

    The triangle is replaced by its three corner copies of half its side, so generation `n` has `3^n` triangles of
    total area `(3/4)^n` times the first, and the gasket has the similarity dimension `log 3 / log 2`, about 1.585. The
    square is replaced by the eight copies of a third of its side around its middle, so generation `n` has `8^n`
    squares of area `(8/9)^n`, and the carpet has the dimension `log 8 / log 3`, about 1.893.

    Args:
        shape:
            `'triangle'` or `'carpet'`.

        filled:
            Whether the geometry visual draws the polygons as a filled image instead of outlines.

    Returns:
        system:
            The system, with the patch `sierpinski` and its number of polygons `sierpinski.size`.

    Raises:
        ValidationError:
            The shape is unknown.

    Examples:
        >>> from skcomplex.systems import create_simulation
        >>> simulation = create_simulation('sierpinski_substitution').prepare().run(until=2)
        >>> simulation.trajectory_.values('sierpinski.size').tolist()
        [1.0, 3.0, 9.0]
    """
    prototile, _ = SierpinskiElement._sierpinski_tiles(shape)
    return System(elements=[SierpinskiElement('sierpinski', shape=shape, prototiles=[prototile], filled=filled)])