---
slug: 2024/11/playing-with-fire
title: "Playing with fire: The fractal flame algorithm"
date: 2024-11-15 12:00:00
authors: [bspeice]
tags: []
---


Wikipedia describes fractal flames fractal flames as:

> a member of the iterated function system class of fractals

It's a bit tedious, but technically correct. I choose to think of them a different way: beauty in mathematics.

import isDarkMode from '@site/src/isDarkMode'
import banner from '../banner.png'

<center>
    <img src={banner} style={{filter: isDarkMode() ? '' : 'invert(1)'}}/>
</center>

<!-- truncate -->

I don't remember exactly when I first learned about fractal flames, but I do remember becoming entranced by the images they created.
I also remember their unique appeal to my young engineering mind; this was an art form I could participate in.

The original [Fractal Flame Algorithm paper](https://flam3.com/flame_draves.pdf) describing their structure was too much
for me to handle at the time (I was ~12 years old), so I was content to play around and enjoy the pictures.
But the desire to understand it stuck around. Now, with a graduate degree under my belt, maybe I can make some progress.

This guide is my attempt to explain fractal flames in a way that younger me &mdash; and others interested in the art &mdash;
can understand without too much prior knowledge.

---

## Iterated function systems

:::note
This post covers section 2 of the Fractal Flame Algorithm paper
:::

As mentioned, fractal flames are a type of "[iterated function system](https://en.wikipedia.org/wiki/Iterated_function_system),"
or IFS. The formula for an IFS is short, but will take some time to unpack:

$$
S = \bigcup_{i=0}^{n-1} F_i(S)
$$

### Solution set

First, $S$. $S$ is the set of points in two dimensions (in math terms, $S \in \mathbb{R}^2$)
that represent a "solution" of some kind to our equation.
Our goal is to find all points in $S$, plot them, and display that image.

For example, if we say $S = \{(0,0), (1, 1), (2, 2)\}$, there are three points to plot:

import {VictoryChart, VictoryTheme, VictoryScatter, VictoryLegend} from "victory";
export const simpleData = [
    {x: 0, y: 0},
    {x: 1, y: 1},
    {x: 2, y: 2}
]

<VictoryChart theme={VictoryTheme.clean}>
    <VictoryScatter data={simpleData} size={5} style={{data: {fill: "blue"}}}/>
</VictoryChart>

However, this is a pretty boring image. With fractal flames, rather than listing individual points,
we use functions to describe which points are part of the solution. This means there are an infinite
number of points, but if we find _enough_ points to plot, we'll end up with a nice picture.
And if we change the functions, our solution changes, and we'll get a new picture.

### Transform functions

Second, the $F_i(S)$ functions, also known as "transforms."
At their most basic, each $F_i$ takes in a 2-dimensional point and gives back a new point
(in math terms, $F_i \in \mathbb{R}^2 \rightarrow \mathbb{R}^2$).
While you could theoretically use any function, we'll focus on a specific kind of function
known as an "[affine transformation](https://en.wikipedia.org/wiki/Affine_transformation)."

The general form of an affine transformation is:

$$
F_i(a_i \cdot x + b_i \cdot y + c_i, d_i \cdot x + e_i \cdot y + f_i)
$$

import transformSource from "!!raw-loader!../src/transform"
import CodeBlock from '@theme/CodeBlock'

<CodeBlock language="typescript">{transformSource}</CodeBlock>

The parameters ($a_i$, $b_i$, etc.) are values we get to choose.
For example, we can represent a "shift" function like this:

$$
\begin{align*}
a &= 1 \\
b &= 0 \\
c &= 0.5 \\
d &= 0 \\
e &= 1 \\
f &= 1.5 \\
F_{shift}(x,y) &= (1 \cdot x + 0 \cdot y + 0.5, 0 \cdot x + 1 \cdot y + 0.5) \\
F_{shift}(x, y) &= (x + 0.5, y + 0.5)
\end{align*}
$$

Applying this function to our original points will give us a new set of points:

import {applyCoefs} from "../src/transform"

export const coefs = {a: 1, b: 0, c: 0.5, d: 0, e: 1, f: 1.5}
export const toData = ([x, y]) => ({x, y})

export const shiftData = simpleData.map(({x, y}) => toData(applyCoefs(x, y, coefs)))

<VictoryChart theme={VictoryTheme.clean}>
    <VictoryScatter data={simpleData} size={5} style={{data: {fill: "blue"}}}/>
    <VictoryScatter data={shiftData} size={5} style={{data: {fill: "orange"}}}/>
    <VictoryLegend
        data={[
            {name: "(x,y)", symbol: {fill: "blue"}},
            {name: "F(x,y)", symbol: {fill: "orange"}}
        ]}
        orientation={"vertical"}
        x={75}
        y={10}
    />
</VictoryChart>

Fractal flames use more complex functions, but they all start with this structure.

### Fixed set

With those definitions in place, let's revisit the initial problem:

$$
S = \bigcup_{i=0}^{n-1} F_i(S)
$$

Or, to put it in English, we would get something like this:

> Our solution, $S$, is the union of all sets produced by applying each function, $F_i$,
> to points in the solution.

There's just one small problem: to find the solution, we must apply these functions to points
we know are in the solution. But how do we know which points are in the solution to start with?

John E. Hutchinson provides an answer in the [original paper](https://maths-people.anu.edu.au/~john/Assets/Research%20Papers/fractals_self-similarity.pdf)
explaining the mathematics of iterated function systems:

> Furthermore, $S$ is compact and is the closure of the set of fixed points $s_{i_1...i_p}$
> of finite compositions $F_{i_1...i_p}$ of members of $F$.

:::note
I've tweaked the wording slightly to match the conventions in the Fractal Flame paper
:::

Before your eyes glaze over, let's unpack this explanation:

- **$S$ is [compact](https://en.wikipedia.org/wiki/Compact_space)...**: All points in our solution will be in a finite range
- **...and is the [closure](https://en.wikipedia.org/wiki/Closure_(mathematics)) of the set of [fixed points](https://en.wikipedia.org/wiki/Fixed_point_(mathematics))**:
  Applying our functions to these points does not change them
- **...of finite compositions $F_{i_1...i_p}$ of members of $F$**: By composing our functions (that is,
  using the output of one function as input to the next function), we will arrive at the points we care about

Thus, by applying the functions in our system to "fixed points," we will find the other points we care about.

However, this is all a bit vague, so let's work through an example.

<details>
    <summary>If you want a bit more math first...</summary>

    ...then there are some details worth mentioning that I've glossed over so far.

    First, the Hutchinson paper requires that the functions $F_i$ be _contractive_ tor the solution set to exist.
    That is, applying the function to a point must bring it closer to other points. However, as the Fractal Flame
    algorithm demonstrates, we only need functions to be contractive _on average_. At worst, the system will
    degenerate and produce a bad image.

    Second, we're focused $\mathbb{R}^2$ because we're generating images, but the Hutchinson paper
    allows for arbitrary dimensions - which means you could also have 3-dimensional fractal flames.

    TODO: Mention attractors? https://en.wikipedia.org/wiki/Chaos_game
</details>

## Sierpinski's gasket

The Fractal Flame paper gives us three functions we can use for our function system:

$$
F_0(x, y) = \left({x \over 2}, {y \over 2} \right) \\
~\\
F_1(x, y) = \left({{x + 1} \over 2}, {y \over 2} \right) \\
~\\
F_2(x, y) = \left({x \over 2}, {{y + 1} \over 2} \right)
$$

### The chaos game

Next, how do we find the "fixed points" we mentioned earlier? The paper lays out an algorithm called the "[chaos game](https://en.wikipedia.org/wiki/Chaos_game)"
that will give us points in the solution set.

$$
\begin{align*}
&(x, y) = \text{random point in the bi-unit square} \\
&\text{iterate } \{ \\
&\hspace{1cm} i = \text{random integer from 0 to } n - 1 \\
&\hspace{1cm} (x,y) = F_i(x,y) \\
&\hspace{1cm} \text{plot}(x,y) \text{ if iterations} > 20 \\
\}
\end{align*}
$$

:::note
In effect, the chaos game algorithm implements the "finite compositions of $F_{i_1..i_p}$ mentioned earlier.
:::

Now, let's turn this into code, one piece at a time.

First, the "bi-unit square" is the range $[-1, 1]$. We can :

import biunitSource from '!!raw-loader!../src/randomBiUnit'

<CodeBlock language="typescript">{biunitSource}</CodeBlock>

Next, we need to choose a random integer from $0$ to $n - 1$:

import randintSource from '!!raw-loader!../src/randomInteger'

<CodeBlock language="typescript">{randintSource}</CodeBlock>

### Plotting

Finally, implementing the `plot` function. This blog series
is designed to be interactive, so everything shows
real-time directly in the browser. As an alternative,
software like `flam3` an Apophysis can also save an image.

To display the results, we'll use the [Canvas API](https://developer.mozilla.org/en-US/docs/Web/API/Canvas_API).
This allows us to manipulate individual pixels an image,
and display it on screen.

First, we need to convert from Fractal Flame coordinates to pixel coordinates.
To simplify things, we'll assume that we're plotting a square image
with range $[0, 1]$ for both $x$ and $y$:

import cameraSource from "!!raw-loader!./cameraGasket"

<CodeBlock language="typescript">{cameraSource}</CodeBlock>

Next, we'll use an [`ImageData` object](https://developer.mozilla.org/en-US/docs/Web/API/ImageData)
to store the pixel data.
Each pixel in the image on screen has a corresponding index in the `data` array.
To plot our image, we set that pixel to be black:

import plotSource from '!!raw-loader!./plot'

<CodeBlock language="typescript">{plotSource}</CodeBlock>

Putting it all together, we have our first image:

import Playground from '@theme/Playground'
import Scope from './scope'

import chaosGameSource from '!!raw-loader!./chaosGame'

<Playground scope={Scope} noInline={true}>{chaosGameSource}</Playground>

<hr/>

<small>
The image here is slightly different than the one in the paper.
I think the paper has an error, so I'm choosing to plot the image
    in a way that's consistent with [`flam3` itself](https://github.com/scottdraves/flam3/blob/7fb50c82e90e051f00efcc3123d0e06de26594b2/rect.c#L440-L441).
</small>

### Weights

Finally, we'll introduce a "weight" ($w_i$) for each function that controls how often we choose
that function during the chaos game relative to each other function.

For Sierpinski's Gasket, we start with equal weighting,
but you can see how changing the function weights affects the image below:

import randomChoiceSource from '!!raw-loader!../src/randomChoice'

<CodeBlock language={'typescript'}>{randomChoiceSource}</CodeBlock>

import chaosGameWeightedSource from "!!raw-loader!./chaosGameWeighted";

<CodeBlock language={'typescript'}>{chaosGameWeightedSource}</CodeBlock>

import GasketWeighted from "./GasketWeighted";
import {SquareCanvas} from "../src/Canvas";

<SquareCanvas><GasketWeighted/></SquareCanvas>

## Summary

Studying the foundations of fractal flames is challenging,
but we now have an understanding of both the mathematics
and implementation of iterated function systems.

In the next post, we'll study the first innovation that fractal flames
bring: variations.