バーンズリーのシダ

**バーンズリーのシダ（Barnsley fern）**は、イギリスの数学者にちなんで名付けられたフラクタルで、その著書『*Fractals Everywhere』*で初めて記述された。チャセンシダの1種であるミヤマカワリシダ（）を模したものである。

## 歴史
シダは、自己相似図形の代表例の1つである。すなわち、数学的に生成されたパターンであり、そのパターンは任意に拡大または縮小しても同一の形状が現れる自己相似性を示す。シェルピンスキーの三角形のように、バーンズリーのシダは、コンピューターで数式を繰り返し使用することで、視覚的に美しい図形を構築する方法を示している。バーンズリーの1988年の本『*Fractals Everywhere*』は、 ジョージア工科大学数学部の*フラクタル幾何学*と呼ばれるコースにおける学部生と大学院生に対する講義内容に基づいている。この本の出版後、*フラクタル測定理論*と呼ばれる第二のコースが設置された。バーンズリーの作品は、数学モデルで自然を模倣しようとするグラフィックアーティストのインスピレーションの源となっている。

バーンズリーによって開発されたシダのコードは、フラクタルを作成するための反復関数系（IFS：Iterated Function System）の例である。これはに従う。彼はフラクタルを使用して、科学技術におけるさまざまな現象、とりわけ植物構造、についてモデル化している。 —マイケル・バーンズリー*ら* 

## シダの生成

バーンズリーのシダは、4つのアフィン変換を使用する。変換の式は次のとおりである。

: 
f(x,y) = \begin{bmatrix}a & b \\ c & d \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} + \begin{bmatrix} e \\ f \end{bmatrix}

あるいは
: 
f(x,y) = \begin{bmatrix}a & b & e \\ c & d & f \\ 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} x \\ y \\ 1 \end{bmatrix}


バーンズリーは、前述した著書で、シダを生成する*IFS*コードを、以下の表のように、値の行列の形で示した。下表では、列「a」から「f」は前記の等式f(x,y)の係数を示しており、「p」はそれぞれの値の組に対する確率因子を表している。

上記の表は、次の4つの変換に対応している。

: 
f_1(x,y) = \begin{bmatrix} \ 0.00 & \ 0.00 \ \\ 0.00 & \ 0.16 \end{bmatrix} \begin{bmatrix} \ x \\ y \end{bmatrix}


: 
f_2(x,y) = \begin{bmatrix} \ 0.85 & \ 0.04 \ \\ -0.04 & \ 0.85 \end{bmatrix} \begin{bmatrix} \ x \\ y \end{bmatrix} + \begin{bmatrix} \ 0.00 \\ 1.60 \end{bmatrix}


: 
f_3(x,y) = \begin{bmatrix} \ 0.20 & \ -0.26 \ \\ 0.23 & \ 0.22 \end{bmatrix} \begin{bmatrix} \ x \\ y \end{bmatrix} + \begin{bmatrix} \ 0.00 \\ 1.60 \end{bmatrix}


: 
f_4(x,y) = \begin{bmatrix} \ -0.15 & \ 0.28 \ \\ 0.26 & \ 0.24 \end{bmatrix} \begin{bmatrix} \ x \\ y \end{bmatrix} + \begin{bmatrix} \ 0.00 \\ 0.44 \end{bmatrix}


### コンピューターによる生成

バーンズリーのシダは、理論的にはペンとグラフ用紙を使用して手でプロットできるが、必要な反復回数は数万回に達するので、実際にはコンピューターを使用する必要がある。バーンズリーのシダの多くのさまざまな計算モデルは、現代の数学者に人気がある。バーンズリーの係数行列を使用し、数式が正しくプログラムされている限り、同じ形状をもったシダが生成される。

最初に描画される点は原点（ *x* 0 = 0、 *y* 0 = 0）にあり、次の4つの座標変換ƒ1, ƒ2, ƒ3, ƒ4のいずれかを設定した発生確率pでランダムに適用することにより、新しい点を反復的に計算する 。すなわち、どの座標変換を使用するか確率に応じた乱数範囲を設定しておき、コンピューターで乱数を発生させ、乱数の値に応じてどの座標変換を使用するのか決め、変換計算を反復的に実行する。これは、1988年にバーンズリーがランダム・アルゴリズムと名付けた初めての例であり、現在ではカオスゲーム、確率論的反復関数系とも呼ばれている。

*ƒ*1 

: *x* *n* + 1 = 0 

: *y* *n* + 1 = 0.16 *y* *n* 。

座標変換ƒ1は1％の確率(p=0.01)で実行される。この変換は、任意の点を、シダの茎基部の最初の線分にマッピングするだけである。フラクタル図のこの部分は、反復の過程で最初に完了する部分である。 

*ƒ*2

: *x* *n* + 1 = 0.85 *x* *n* + 0.04 *y* *n* 

: *y* *n* + 1 = −0.04 *x* *n* + 0.85 *y* *n* + 1.6。

この座標変換ƒ2は、85％の確率(p=0.85)で実行され、図の赤い三角形で表される葉の中の任意のポイントを、青い三角形で表される反対側の小さな葉の中の点にマッピングする。

*ƒ*3

: *x* *n* + 1 = 0.2 *x* *n* &nbsp; − &nbsp; 0.26 &nbsp; *y* *n* 

: *y* *n* + 1 = 0.23 *x* *n* &nbsp; + &nbsp; 0.22 *y* *n* &nbsp; + &nbsp; 1.6。

この座標変換ƒ3は、7％の確率(p=0.07)で実行され、青い三角形で表される葉（または羽片）の中の任意の点を、茎を交差し交互に対応する三角形内にマッピングする（反転する）。

*ƒ*4

: *x* *n* + 1 = −0.15 *x* *n* + 0.28 *y* *n* 

: *y* *n* + 1 = 0.26 *x* *n* + 0.24 *y* *n* + 0.44。

この座標変換ƒ4は、7％の確率(p=0.07)で実行され、青い三角形で表される葉（または羽片）の中の任意の点を、茎を交差し交互に対応する三角形内にマッピングする（反転せずに）。

最初の座標変換は茎を描画する。2番目は茎と下葉を連続コピーして、完全なシダを生成する。3番目は左下葉を描き、4番目は、右下の葉を描く。IFSの再帰的な性質により、全体が各葉の部分の大きな複製となっている。シダ全体が-2.1820 
import turtle
import random

pen = turtle.Turtle()
pen.speed(20)
pen.color('green')
pen.penup()

x = 0
y = 0
for n in range(100000):
    pen.goto(85*x, 57*y-275)  # 57 is to scale the fern and -275 is to start the drawing from the bottom.
    pen.pendown()
    pen.dot()
    pen.penup()
    r = random.random()  # to get probability
    r = r * 100
    xn = x
    yn = y
    if r 

### R

```r
# Barnsley's Fern

# create function of the probability and the current point
fractal_fern2 <- function(x, p){
  if (p <= 0.01) {
    m <- matrix(c(0, 0, 0, .16), 2, 2)
    f <- c(0, 0)
  } else if (p <= 0.86) {
    m <- matrix(c(.85, -.04, .04, .85), 2, 2)
    f <- c(0, 1.6)
  } else if (p <= 0.93) {
    m <- matrix(c(.2, .23, -.26, .22), 2, 2)
    f <- c(0, 1.6)
  } else {
    m <- matrix(c(-.15, .26, .28, .24), 2, 2)
    f <- c(0, .44)
  }
  m %*% x + f
}

# how many reps determines how detailed the fern will be
reps <- 10000

# create a vector with probability values, and a matrix to store coordinates
p <- runif(reps)

# initialise a point at the origin
coords <- c(0, 0)

# compute Fractal Coordinates
m <- Reduce(fractal_fern2, p, accumulate = T, init = coords)
m <- t(do.call(cbind, m))

# Create plot
plot(m, type = "p", cex = 0.1, col = "darkgreen",
     xlim = c(-3, 3), ylim = c(0, 10), 
     xlab = NA, ylab = NA, axes = FALSE)
```


### Processing

```java
/*
    Barnsley Fern for Processing 3.4
*/

// declaring variables x and y
float x, y;

// creating canvas
void setup() {
  size(600, 600);
  background(255);
}

/* setting stroke,  mapping canvas and then
   plotting the points */
void drawPoint() {
  stroke(34, 139, 34);
  strokeWeight(1);
  float px = map(x, -2.1820, 2.6558, 0, width);
  float py = map(y, 0, 9.9983, height, 0);
  point(px, py);
}

/* algorithm for calculating value of (n+1)th
   term of x and y based on the transformation
   matrices */
void nextPoint() {
  float nextX, nextY;
  float r = random(1);
  if (r < 0.01) {
    nextX =  0;
    nextY =  0.16 * y;
  } else if (r < 0.86) {
    nextX =  0.85 * x + 0.04 * y;
    nextY = -0.04 * x + 0.85 * y + 1.6;
  } else if (r < 0.93) {
    nextX =  0.20 * x - 0.26 * y;
    nextY =  0.23 * x + 0.22 * y + 1.6;
  } else {
    nextX = -0.15 * x + 0.28 * y;
    nextY =  0.26 * x + 0.24 * y + 0.44;
  }
  x = nextX;
  y = nextY;
}

/* iterate the plotting and calculation
   functions over a loop */
void draw() {
  for (int i = 0; i < 100; i++) {
    drawPoint();
    nextPoint();
  }
}
```


### P5.JS

```javascript
let x = 0;
let y = 0;

function setup() {
 createCanvas(600, 600);
 background(0);
}

//range −2.1820 < x < 2.6558 and 0 ≤ y < 9.9983.
function drawPoint() {
 stroke(255);
 strokeWeight(1);
 let px = map(x, -2.1820, 2.6558, 0, width);
 let py = map(y, 0, 9.9983, height, 0);
 point(px, py);
}

function nextPoint() {
 let nextX;
 let nextY;

 let r = random(1);

 if (r < 0.01) {
  //1
  nextX = 0;
  nextY = 0.16 * y;
 } else if (r < 0.86) {
  //2
  nextX = 0.85 * x + 0.04 * y;
  nextY = -0.04 * x + 0.85 * y + 1.60;
 } else if (r < 0.93) {
  //3
  nextX = 0.20 * x + -0.26 * y;
  nextY = 0.23 * x + 0.22 * y + 1.60;
 } else {
  //4
  nextX = -0.15 * x + 0.28 * y;
  nextY = 0.26 * x + 0.24 * y + 0.44;
 }

 x = nextX;
 y = nextY;
}

function draw() {
 for (let i = 0; i < 1000; i++) {
  drawPoint();
  nextPoint();
 }
}
```


### JavaScript (HTML5)

```html
<canvas id="canvas" height="700" width="700">
</canvas>

<script>
    let canvas;
    let canvasContext;

    let move;

    let fps = 250;

    let pointColor = "green";
    let backgroundColor = "black";
    let pointRadius = 1;

    let x = 0;
    let y = 0;

    window.onload = function () {
        canvas = document.getElementById("canvas");
        canvasContext = canvas.getContext('2d');

        canvasContext.fillStyle = backgroundColor;
        canvasContext.fillRect(0, 0, canvas.width, canvas.height);

        move = setInterval(doIt, 1000/fps);
    };

    let doIt = function () {
        for (let i = 0; i < 20; i++)
            moveAll()
    };

    function moveAll() {

        let nextX, nextY;
        let r = Math.random();
        if (r < 0.01) {
            nextX =  0;
            nextY =  0.16 * y;
        } else if (r < 0.86) {
            nextX =  0.85 * x + 0.04 * y;
            nextY = -0.04 * x + 0.85 * y + 1.6;
        } else if (r < 0.93) {
            nextX =  0.20 * x - 0.26 * y;
            nextY =  0.23 * x + 0.22 * y + 1.6;
        } else {
            nextX = -0.15 * x + 0.28 * y;
            nextY =  0.26 * x + 0.24 * y + 0.44;
        }

        // Scaling and positioning
        let plotX = canvas.width * (x + 3) / 6;
        let plotY = canvas.height - canvas.height * ((y + 2) / 14);

        drawFilledCircle(plotX, plotY, pointRadius, pointColor);

        x = nextX;
        y = nextY;

    }
    const drawFilledCircle = (centerX, centerY, radius, color) => {
        canvasContext.beginPath();
        canvasContext.fillStyle = color;
        canvasContext.arc(centerX, centerY, radius, 0, 2 * Math.PI, true);
        canvasContext.fill();
    };
</script>
```


### QBasic

```qbasic
SCREEN 12
WINDOW (-5, 0)-(5, 10)
RANDOMIZE TIMER
COLOR 10
DO
    SELECT CASE RND
        CASE IS < .01
            nextX = 0
            nextY = .16 * y
        CASE .01 TO .08
            nextX = .2 * x - .26 * y
            nextY = .23 * x + .22 * y + 1.6
        CASE .08 TO .15
            nextX = -.15 * x + .28 * y
            nextY = .26 * x + .24 * y + .44
        CASE ELSE
            nextX = .85 * x + .04 * y
            nextY = -.04 * x + .85 * y + 1.6
    END SELECT
    x = nextX
    y = nextY
    PSET (x, y)
LOOP UNTIL INKEY$ = CHR$(27)
```


### VBA（CorelDraw）

```basic
Sub Barnsley()

Dim iEnd As Long
Dim i As Long
Dim x As Double
Dim y As Double
Dim nextX As Double
Dim nextY As Double
Dim sShapeArray() As Shape
Dim dSize As Double
Dim sColor As String


    dSize = 0.01        ' Size of the dots
    sColor = "0,0,100"  ' RGB color of dots, value range 0 to 255
    iEnd = 5000       ' Number of iterations
    ReDim sShapeArray(iEnd)
                        ' In Corel, each object drawn requires a variable name of its own

    Randomize           ' Initialize the Rnd function

    For i = 0 To iEnd   ' Iterate ...
        Select Case Rnd
            Case Is < 0.01
                        ' f1 = Draw stem
                nextX = 0
                nextY = 0.16 * y
            Case 0.01 To 0.08
                        ' f3
                nextX = 0.2 * x - 0.26 * y
                nextY = 0.23 * x + 0.22 * y + 1.6
            Case 0.08 To 0.15
                        ' f4
                nextX = -0.15 * x + 0.28 * y
                nextY = 0.26 * x + 0.24 * y + 0.44
            Case Else
                        ' f2
                nextX = 0.85 * x + 0.04 * y
                nextY = -0.04 * x + 0.85 * y + 1.6
        End Select
        x = nextX
        y = nextY
        Set sShapeArray(i) = ActiveLayer.CreateEllipse2(x + 2.5, y + 0.5, dSize)
        sShapeArray(i).Style.StringAssign "{""fill"":{""primaryColor"":""RGB255,USER," & sColor & ",100,00000000-0000-0000-0000-000000000000"",""secondaryColor"":""RGB255,USER,255,255,255,100,00000000-0000-0000-0000-000000000000"",""type"":""1"",""fillName"":null},""outline"":{""width"":""0"",""color"":""RGB255,USER,0,0,0,100,00000000-0000-0000-0000-000000000000""},""transparency"":{}}"
        DoEvents
    Next

 End Sub
```


### Amola

```abap
addpackage("Forms.dll")

set("x", 0)
set("y", 0)
set("width", 600)
set("height", 600)

method setup()
	createCanvas(width, height)
	rect(0, 0, 600, 600, color(0, 0, 0))
end

method drawPoint()
    set("curX", div(mult(width, add(x, 3)), 6))
    set("curY", sub(height, mult(height, div(add(y, 2), 14))))
    set("size", 1)
	//log(curX)
	//log(curY)
	rect(round(curX - size / 2), round(curY - size / 2), round(curX + size / 2), round(curY + size / 2), color(34, 139, 34))
end

method nextPoint()
	set("nextX", 0)
	set("nextY", 0)
	set("random", random(0, 100))
	if(random < 1)
		set("nextX", 0)
		set("nextY", 0.16 * y)
	end
	else
		if(random < 86)
			set("nextX", 0.85 * x + 0.04 * y)
			set("nextY", -0.04 * x + 0.85 * y + 1.6)
		end
		else
			if(random < 93)
				set("nextX", 0.2 * x - 0.26 * y)
				set("nextY", 0.23 * x + 0.22 * y + 1.6)
			end
			else
				set("nextX", -0.15 * x + 0.28 * y)
				set("nextY", 0.26 * x + 0.24 * y + 0.44)
			end
		end
	end

	set("x", nextX)
	set("y", nextY)
end

setup()
while(true)
	drawPoint()
	nextPoint()
end
```


### 表計算ソフトによる計算例
プログラム言語の構文と同じ計算を、以下のように表計算ソフトの関数を使用することでも実現できる。

なお、B8,C8,D8のセルには以下のような複数条件判定の関数(ネスティング参照)を入力する。
*B8=IF(A8<($H$2),1,IF(A8<($H$2+$H$3),2,IF(A8<($H$2+$H$3+$H$4),3,4)))
*C8=IF(B8=1,$B$2*C7+$C$2*D7+$F$2,IF(B8=2,$B$3*C7+$C$3*D7+$F$3,IF(B8=3,$B$4*C7+$C$4*D7+$F$4,$B$5*C7+$C$5*D7+$F$5)))
*D8=IF(B8=1,$D$2*C7+$E$2*D7+$G$2,IF(B8=2,$D$3*C7+$E$3*D7+$G$3,IF(B8=3,$D$4*C7+$E$4*D7+$G$4,$D$5*C7+$E$5*D7+$G$5)))
最終8行目をオートフィルで適当な行数だけコピーし、XY散布図とするとバーンズリーのシダのフラクタル図が得られる。各変換式ƒの係数a,b,c,d,e,fと確率pは任意に変更可能である。

各列は以下のような計算を行っている。
*A列	乱数を発生させる。
*B列	乱数をもとに確率pに応じた条件判定を行い、用いる変換ƒを決める。
*C列	先に決めた変換ƒに対応する計算をおこない、Xを求める。
*D列	先に決めた変換ƒに対応する計算をおこない、Yを求める。
*新たなXとYは前の行のXとYの値を使用し、反復的に計算を進める。

## 脚注