Шёл уже последний час этого воскресенья, я уже думал идти спать, но добрый sourcerer прислал мне картинку с моего заброшенного сайта, которую можно увидеть ниже, и текст «красиво!». Эти картинки я рисовал лет пять назад, с помощью т. н. алгоритма времени убегания, но для применимости данного алгоритма, нужно уметь для заданного набора преобразований разбивать плоскость на регионы, тогда я не придумал, как это сделать, и больше к этому алгоритму не возвращался. Но сейчас я сразу сообразил, что делать, и написал Диме: «Сначала Random IFS, потом kNN, а затем Escape-Time Algorithm!»
![](//habrastorage.org/r/w1560/files/553/eef/519/553eef51989c44aaac42f949378c6d0c.png)
Под рукой у меня был только старый нетбук, который мне дали друзья на время, пока мой ноутбук в ремонте. Дима мне ещё что-то говорил, я ему что-то отвечал, но у меня уже в голове писался код, и я искал на нетбуке хоть какой-нибудь компилятор или интерпретатор и нашёл C++ Builder 6! После этого я понял, что утро я встречу наедине с борландовским компилятором. Через пять часов я отправил Диме новых картинок, но он, как нормальный человек, давно спал…
![](//habrastorage.org/files/536/9e7/fa3/5369e7fa36994d92a59b82d19cc83f78.PNG)
Итак, немножко формул. Представим, что у нас есть конечный набор преобразований плоскости Ti: R2 → R2, i = 1, ..., k. Для произвольного множества E определим T(E) = T1(E) ∪… ∪ Tk(E), т. е. подействуем каждым преобразованием на множество E, а результаты объединим. Можно доказать, что если отображения Ti были сжимающими, то последовательность E, T(E), T(T(E)),… сойдётся к некоторому множеству F для любого непустого компактного E. Данная конструкция и известна как система итерируемых функций.
Например, если в качестве непустого компакта взять смайлик, и рассмотреть три преобразования, каждое из которых является композицией сжатия и сдвига в i-ую вершину правильного треугольника, то первые итерации будут выглядеть так, а в пределе получится треугольник Серпинского:
![](https://habrastorage.org/files/c9d/7ef/1f0/c9d7ef1f0c56419ba10e9fe3b4bb173d.PNG)
Обычно вместо прямого вычисления последовательности E, T(E), T(T(E)),… для построения фракталов используют т. н. «игру в хаос», которая заключается в следующем. Выберем произвольную точку z0 на плоскости, далее выберем случайно преобразование Ti1 и вычислим z1= Ti1(z0), далее снова случайно выберем Ti2 и вычислим z1= Ti2(z0), и т. д. Можно показать, что всё будет хорошо, и множество полученных точек будет в некотором смысле приближать множество F, определённое выше. На этот алгоритм я ниже буду ссылаться как на Random IFS.
![](https://habrastorage.org/files/310/080/30f/31008030fffd4cd199ad66d38f96f297.PNG)
Теперь самое время перейти к описанию алгоритма времени убегания. Пусть у нас для начала есть одно преобразование плоскости f. Для каждой точки z плоскости начнём вычислять последовательность f(z), f(f(z)), f(f(f(z))),… до тех пор пока либо число итераций не превысит некоторого заданного числа N, либо пока норма числа z не станет больше некоторого числа B. После этого цвет точки выбираем в соответствии с количеством произведённых итераций.
Если на время представить что наша плоскость является комплексной, а преобразование f(z) равно z2 + c, то в результате работы этого алгоритма мы получим фрактальное множество Жюлиа. Более подробно про это можно прочитать в хабрастатье «Построение множества Жюлиа» хабрапользователя mephistopheies.
![](//habrastorage.org/r/w780q1/files/33f/f7c/0e1/33ff7c0e135e4c478e967d8981742892.jpg)
Пусть теперь у нас есть система итерируемых функций, заданная набором обратимых сжимающих преобразований плоскости T1, ..., Tk. Пусть F — это аттрактор этой системы.
Дополнительно предположим, что множество F можно разбить так, что Ti(F) ∩ Tj(F) = ∅, i != j (это предположение далеко не всегда выполняется). Разобъём всю плоскость R2 на куски A1, ...., Ak так, что Ti(F) является подмножеством Ai для всех i. Теперь определим функцию f(z) кусочно: на множестве Ai положим f(z) = Tk−1(z) для всех i.
Например, для треугольника Серпинского рассмотрим такое разбиение (тут есть небольшие проблемы с тремя точками, но закроем на них глаза).
![](//habrastorage.org/files/48b/d61/5f5/48bd615f52244890a5d8693213ad084d.PNG)
А теперь самый главный вопрос, что получится, если алгоритм времени убегания применить к построенной таким образом функции f?
Давайте посмотрим:
![](//habrastorage.org/r/w1560/files/8c4/ec9/a65/8c4ec9a6550446c6908f185e9c581baf.png)
Получился симпатичный такой треугольник Серпинского!
Оказывается это не случайность. Ещё пару примеров:
![](//habrastorage.org/r/w1560/files/303/0cd/b74/3030cdb7488447c891f6ad73efeb442b.png)
![](//habrastorage.org/r/w1560/files/553/eef/519/553eef51989c44aaac42f949378c6d0c.png)
В этих примерах соответствующее разбиение плоскости не сложно задать аналитически с помощью булевых комбинаций кругов и полуплоскостей, используя метод пристального вглядывания. Но часто простых условий угадать не удаётся. Поэтому вместо угадывания мы научим компьютер определять разбиение самостоятельно. В этом нам поможет метод ближайшего соседа.
А именно, сначала с помощью Random IFS генерируем несколько тысяч точек, при этом для каждой точки запоминаем номер преобразования, с помощью которого она была получена. Затем во время работы EscapeTimeAlgorithm для каждого пикселя определяем область, в которую он попадаем с помощью 1NN.
Например, для такой звёздочки 1NN даёт следующее разбиение на четыре куска:
![](//habrastorage.org/files/728/953/46f/72895346fb294254a5aba2517a3d120e.JPG)
Собирая вместе, получим:
![](//habrastorage.org/files/f64/d5c/8f4/f64d5c8f4c114420ab4490c26000e1f5.PNG)
Ещё несколько картинок.
![](//habrastorage.org/files/cc6/c52/d22/cc6c52d22e774f6286a5172874a1a0b3.PNG)
![](//habrastorage.org/files/536/9e7/fa3/5369e7fa36994d92a59b82d19cc83f78.PNG)
![](//habrastorage.org/files/c7e/7d3/888/c7e7d38888144c51b6467b5ed019a209.PNG)
Вот и всё. Напоследок два замечания.
Во-первых, внимательный читатель возможно задался вопросом, раз фракталы, которые строятся с помощью Random IFS, можно построить с помощью алгоритма времени убегания, то можно ли множество Жюлиа построить с помощью Random IFS? Оказывается можно, нужно просто обратить отображение f(z) = z2 + c, вспомнив, как извлекается корень из комплексного числа. (Правда при применении этого метода для построения изображений множества Жюлиа возникают большие трудности.)
![](//habrastorage.org/files/cec/fb8/439/cecfb843996241f48ea9073681ed655e.PNG)
Во-вторых, в статье было рассказано что происходит, если вы хотите узнать почему так происходит, то рекомендую книгу M. Barnsley «Fractals Everywere».
(Мгновения исходного кода можно найти на гитхабе.)
![](http://habrastorage.org/files/553/eef/519/553eef51989c44aaac42f949378c6d0c.png)
Под рукой у меня был только старый нетбук, который мне дали друзья на время, пока мой ноутбук в ремонте. Дима мне ещё что-то говорил, я ему что-то отвечал, но у меня уже в голове писался код, и я искал на нетбуке хоть какой-нибудь компилятор или интерпретатор и нашёл C++ Builder 6! После этого я понял, что утро я встречу наедине с борландовским компилятором. Через пять часов я отправил Диме новых картинок, но он, как нормальный человек, давно спал…
Итак, немножко формул. Представим, что у нас есть конечный набор преобразований плоскости Ti: R2 → R2, i = 1, ..., k. Для произвольного множества E определим T(E) = T1(E) ∪… ∪ Tk(E), т. е. подействуем каждым преобразованием на множество E, а результаты объединим. Можно доказать, что если отображения Ti были сжимающими, то последовательность E, T(E), T(T(E)),… сойдётся к некоторому множеству F для любого непустого компактного E. Данная конструкция и известна как система итерируемых функций.
Например, если в качестве непустого компакта взять смайлик, и рассмотреть три преобразования, каждое из которых является композицией сжатия и сдвига в i-ую вершину правильного треугольника, то первые итерации будут выглядеть так, а в пределе получится треугольник Серпинского:
Обычно вместо прямого вычисления последовательности E, T(E), T(T(E)),… для построения фракталов используют т. н. «игру в хаос», которая заключается в следующем. Выберем произвольную точку z0 на плоскости, далее выберем случайно преобразование Ti1 и вычислим z1= Ti1(z0), далее снова случайно выберем Ti2 и вычислим z1= Ti2(z0), и т. д. Можно показать, что всё будет хорошо, и множество полученных точек будет в некотором смысле приближать множество F, определённое выше. На этот алгоритм я ниже буду ссылаться как на Random IFS.
z = (0, 0)
for (i = 0; i < maxIterNum; ++i) {
cl = random([p1, ..., pk]) // pi -- вероятность, с которой выбираем преобразование Ti.
z = T[cl](z)
if (i > skipIterNum) { // Первых несколько итераций могут быть достаточно далеко от аттрактора.
draw(z)
}
}
Теперь самое время перейти к описанию алгоритма времени убегания. Пусть у нас для начала есть одно преобразование плоскости f. Для каждой точки z плоскости начнём вычислять последовательность f(z), f(f(z)), f(f(f(z))),… до тех пор пока либо число итераций не превысит некоторого заданного числа N, либо пока норма числа z не станет больше некоторого числа B. После этого цвет точки выбираем в соответствии с количеством произведённых итераций.
for (y = y1; y < y2; y += dy) {
for (x = x1; x < x2; x += dx) {
z = (x, y);
iter = 0;
while (iter < N && ||z|| < B) {
z = f(z)
iter += 1;
}
drawPixel((x, y), color(iter))
}
}
Если на время представить что наша плоскость является комплексной, а преобразование f(z) равно z2 + c, то в результате работы этого алгоритма мы получим фрактальное множество Жюлиа. Более подробно про это можно прочитать в хабрастатье «Построение множества Жюлиа» хабрапользователя mephistopheies.
![](http://habrastorage.org/files/33f/f7c/0e1/33ff7c0e135e4c478e967d8981742892.jpg)
Пусть теперь у нас есть система итерируемых функций, заданная набором обратимых сжимающих преобразований плоскости T1, ..., Tk. Пусть F — это аттрактор этой системы.
Дополнительно предположим, что множество F можно разбить так, что Ti(F) ∩ Tj(F) = ∅, i != j (это предположение далеко не всегда выполняется). Разобъём всю плоскость R2 на куски A1, ...., Ak так, что Ti(F) является подмножеством Ai для всех i. Теперь определим функцию f(z) кусочно: на множестве Ai положим f(z) = Tk−1(z) для всех i.
Например, для треугольника Серпинского рассмотрим такое разбиение (тут есть небольшие проблемы с тремя точками, но закроем на них глаза).
А теперь самый главный вопрос, что получится, если алгоритм времени убегания применить к построенной таким образом функции f?
Давайте посмотрим:
![](http://habrastorage.org/files/8c4/ec9/a65/8c4ec9a6550446c6908f185e9c581baf.png)
Получился симпатичный такой треугольник Серпинского!
Оказывается это не случайность. Ещё пару примеров:
![](http://habrastorage.org/files/808/d11/43f/808d1143fc2442c8bc78d48c9fe1a381.png)
![](http://habrastorage.org/files/303/0cd/b74/3030cdb7488447c891f6ad73efeb442b.png)
![](http://habrastorage.org/files/c7b/010/f87/c7b010f87438440cae6c33e7c5b6eb61.png)
![](http://habrastorage.org/files/553/eef/519/553eef51989c44aaac42f949378c6d0c.png)
В этих примерах соответствующее разбиение плоскости не сложно задать аналитически с помощью булевых комбинаций кругов и полуплоскостей, используя метод пристального вглядывания. Но часто простых условий угадать не удаётся. Поэтому вместо угадывания мы научим компьютер определять разбиение самостоятельно. В этом нам поможет метод ближайшего соседа.
А именно, сначала с помощью Random IFS генерируем несколько тысяч точек, при этом для каждой точки запоминаем номер преобразования, с помощью которого она была получена. Затем во время работы EscapeTimeAlgorithm для каждого пикселя определяем область, в которую он попадаем с помощью 1NN.
Например, для такой звёздочки 1NN даёт следующее разбиение на четыре куска:
![](http://habrastorage.org/files/a8a/2fb/a49/a8a2fba492eb467da348a3fbacbdfa74.png)
Собирая вместе, получим:
points = RandomIFS(Ts)
classifier = kNN(points);
for (y = y1; y < y2; y += dy) {
for (x = x1; x < x2; x += dx) {
z = (x, y)
iter = 0
while (iter < maxIterNum && ||z|| < sqrdBound) {
cl = classifier.getClass(z);
z = T[cl].applyInvert(z);
iter += 1;
}
draw((x, y), color(iter))
}
}
Ещё несколько картинок.
Вот и всё. Напоследок два замечания.
Во-первых, внимательный читатель возможно задался вопросом, раз фракталы, которые строятся с помощью Random IFS, можно построить с помощью алгоритма времени убегания, то можно ли множество Жюлиа построить с помощью Random IFS? Оказывается можно, нужно просто обратить отображение f(z) = z2 + c, вспомнив, как извлекается корень из комплексного числа. (Правда при применении этого метода для построения изображений множества Жюлиа возникают большие трудности.)
x = z0.re
y = z0.im
for (i = 0; i < N; ++i) {
x -= c.re;
y -= c.im;
len = sqrt(x*x + y*y);
ang = atan2(y, x);
if (rand() % 2 == 0) { // Тут нужно что-нибудь и поинтереснее.
x = sqrt(len) * cos(ang / 2);
y = sqrt(len) * sin(ang / 2);
} else {
x = sqrt(len) * cos(ang / 2 + M_PI);
y = sqrt(len) * sin(ang / 2 + M_PI);
}
draw(x, y)
}
Во-вторых, в статье было рассказано что происходит, если вы хотите узнать почему так происходит, то рекомендую книгу M. Barnsley «Fractals Everywere».
(Мгновения исходного кода можно найти на гитхабе.)