IFS(Iterated Function System)反復関数系

# フラクタルの「シダ」

よくフラクタルで登場する「シダ」のグラフィックス。
以下４種のx,yの漸化式が用意されており、確率でどれかを選び、計算を繰り返しながら「点」を描画する。

```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0 \\
y_{n+1}=0.16y_n
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.2x_n-0.26y_n \\
y_{n+1}=0.23x_n+0.22y_n+1.6
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=-0.15x_n+0.28y_n \\
y_{n+1}=0.26y_n+0.24y_n+0.44
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.85x_n+0.04y_n \\
y_{n+1}=-0.04x_n+0.85y_n+1.6
   \end{cases}
\end{equation*}
```


回転行列による計算を
```math
\begin{bmatrix}
x_{n+1} \\
y_{n+1} \\
\end{bmatrix}
=
\begin{bmatrix}
cos\theta & -sin\theta \\
sin\theta & cos\theta \\
\end{bmatrix}
\begin{bmatrix}
x_n \\
y_n \\
\end{bmatrix}
```
とすると
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=x_ncos\theta-y_nsin\theta \\
y_{n+1}=x_nsin\theta+y_ncos\theta
   \end{cases}
\end{equation*}
```
それで、最後の式で回転角度調べてみようと、
```math
arccos\theta=0.85
```
```math
\theta=0.554811033rad=0.00968327925degree
```
```math
0.85^2+0.04^2=0.7241
```
単純な回転じゃなかった。せん断が入ってるのか？

```java
size(256, 512);
float x=0, y=0;
float x1, y1;
translate(width/2, 0);
for (int n=0; n<100000; n++) {
  float r = random(100);
  if (r<1) {
    x1 = 0;
    y1 = 0.16*y;
  } else if (r<8) {
    x1 =  0.2*x-0.26*y;
    y1 = 0.23*x+0.22*y+1.6;
  } else if (r<15) {
    x1 = -0.15*x+0.28*y;
    y1 =  0.26*x+0.24*y+0.44;
  } else {
    x1 =  0.85*x+0.04*y;
    y1 = -0.04*x+0.85*y+1.6;
  }
  x = x1;
  y = y1;
  stroke(r);
  point(x*40, y*40);
}
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/510257/9d3a3cbd-683d-082d-3808-1830ea190568.png)

# シェルピンスキーのギャスケット
0.5倍
0.5倍して(0.25, √3/4)の移動
0.5倍して(0.5, 0)の移動
の3種のx,yの漸化式を用意し、以下同様。
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.5x_n \\
y_{n+1}=0.5y_n
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.5x_n+0.25 \\
y_{n+1}=0.5y_n+\frac{\sqrt3}{4}
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.5x_n+0.5 \\
y_{n+1}=0.5y_n
   \end{cases}
\end{equation*}
```

```java
size(512,512);
float Xn=0, Yn=0;
float X1, Y1;

for (int n=0;n<100000; n++) {
  float r = random(100);
  if (r<33) {
    X1 = 0.5*Xn;
    Y1 = 0.5*Yn;
  } else if (r<66) {
    X1=0.5*Xn+0.25;
    Y1=0.5*Yn+sqrt(3f)/4f;
  } else{
    X1=0.5*Xn+0.5;
    Y1=0.5*Yn;
  }
  Xn=X1;
  Yn=Y1;
  stroke(r);
  point(Xn*512, Yn*512);
}
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/510257/e534097f-a11d-35d2-7fb7-88eb829d1bb6.png)


# コッホ曲線
```java
size(512,512);
float Xn=0, Yn=0;
float X1, Y1;
translate(0, height/2);
for (int n=0;n<10000; n++) {
  float r = random(100);
  float s = 1f/3f;
  float t = radians(60);
  if (r<25) {
    X1 = s*Xn;
    Y1 = s*Yn;
  } else if (r<50) {
    X1 = s*cos(t)*Xn - s*sin(t)*Yn + s;
    Y1 = s*sin(t)*Xn + s*cos(t)*Yn;
  } else if (r<75) {
    X1 = s*cos(t)*Xn + s*sin(t)*Yn + cos(t);
    Y1 =-s*sin(t)*Xn + s*cos(t)*Yn + s*sin(t);
  } else{
    X1 = s*Xn + 2*s;
    Y1 = s*Yn;
  }
  Xn = X1;
  Yn = Y1;
  point(Xn*512, Yn*512);
}
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/510257/13751ff0-9432-b07e-ac39-763f1a239836.png)

# ドラゴン曲線
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.5x_n-0.5y_n \\
y_{n+1}=0.5x_n+0.5y_n
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=-0.5x_n-0.5y_n+1 \\
y_{n+1}=0.5x_n-0.5y_n
   \end{cases}
\end{equation*}
```

```java
size(512,512);
float Xn=0, Yn=0;
float X1, Y1;
translate(width/4, height/2);
for (int n=0;n<200000; n++) {
  float r = random(100);
  if (r<50) {
    X1 = 0.5*Xn-0.5*Yn;
    Y1 = 0.5*Xn+0.5*Yn;
  } else{
    X1 =-0.5*Xn-0.5*Yn + 1;
    Y1 = 0.5*Xn-0.5*Yn;
  }
  Xn = X1;
  Yn = Y1;
  stroke(r);
  point(Xn*256, Yn*256);
}
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/510257/6ccf5f10-b920-89d0-0c9d-b0033db3a82d.png)


# レヴィドラゴン曲線
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.5x_n-0.5y_n \\
y_{n+1}=0.5x_n+0.5y_n
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=-0.5x_n+0.5y_n+1 \\
y_{n+1}=0.5x_n+0.5y_n
   \end{cases}
\end{equation*}
```

```java
size(512,512);
float Xn=0, Yn=0;
float X1, Y1;
translate(width/4, height/2);
for (int n=0;n<100000; n++) {
  float r = random(100);
  if (r<50) {
    X1 = 0.5*Xn-0.5*Yn;
    Y1 = 0.5*Xn+0.5*Yn;
  } else{
    X1 =-0.5*Xn+0.5*Yn + 1;
    Y1 = 0.5*Xn+0.5*Yn;
  }
  Xn = X1;
  Yn = Y1;
  stroke(r);
  point(Xn*256, Yn*256);
}
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/510257/504fe293-6bc2-1b5b-693d-3e934c67453d.png)


# ツインドラゴン
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=0.5x_n-0.5y_n \\
y_{n+1}=0.5x_n+0.5y_n
   \end{cases}
\end{equation*}
```
```math
\begin{equation*}
   \begin{cases}
x_{n+1}=-0.5x_n+0.5y_n+1 \\
y_{n+1}=-0.5x_n-0.5y_n
   \end{cases}
\end{equation*}
```

```java
size(512,512);
float Xn=0, Yn=0;
float X1, Y1;
translate(width/4, height/2);
for (int n=0;n<1000000; n++) {
  float r = random(100);
  if (r<50) {
    X1 = 0.5*Xn-0.5*Yn;
    Y1 = 0.5*Xn+0.5*Yn;
  } else{
    X1 =-0.5*Xn+0.5*Yn + 1;
    Y1 =-0.5*Xn-0.5*Yn;
  }
  Xn = X1;
  Yn = Y1;
  stroke(r);
  point(Xn*256, Yn*256);
}
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/510257/2f2a3da0-af9c-eb27-1280-82461ec760b3.png)



#　参考
反復関数系

https://ja.wikipedia.org/wiki/%E5%8F%8D%E5%BE%A9%E9%96%A2%E6%95%B0%E7%B3%BB

Iterated function system

https://en.wikipedia.org/wiki/Iterated_function_system