Перейти к содержанию

Множество Мандельброта

завершено на 25%
Материал из Викиучебника — открытых книг для открытого мира
Множество Мандельброта
Множество Мандельброта
Логотип Википедии
Логотип Википедии

В Википедии имеется статья по теме «Множество Мандельброта»

Мно́жество Мандельбро́та — это множество таких точек c на комплексной плоскости, для которых рекуррентное соотношение при задаёт ограниченную последовательность. То есть, это множество таких c, для которых существует такое действительное R, что неравенство |zn|<R выполняется при всех натуральных n.


Построение множества

[править]

Несложно доказать, что как только модуль zn окажется больше 2 (или, в терминах действительной и мнимой частей, xn2+yn2>4), последовательность станет стремиться к бесконечности. В случае |c|≤2 это можно доказать с помощью метода математической индукции. При |c|>2 точка c заведомо не принадлежит множеству Мандельброта, что также можно вывести методом индукции, используя равенство z0=0. (Хотя в этом случае может существовать другое z0, для которого соответствующая последовательность ограничена по модулю, но для некоторого n выполняется неравенство |zn|>2.)

Сравнение |zn| с этим числом (в англоязычной литературе его называют «bail-out») позволяет выделять точки, не попадающие внутрь множества. Для точек, лежащих внутри множества, последовательность не будет иметь тенденции к бесконечности и никогда не достигнет этого числа, поэтому после определённого числа итераций расчёт необходимо принудительно завершить. Максимальное число итераций, после которых число считается попавшим внутрь множества, задается в программе.

Изображение, полученное таким способом, является лишь приближением к реальному множеству Мандельброта. Более качественные результаты можно получать, увеличивая максимальное количество итераций, однако при этом пропорционально вырастает и время расчётов.

Примеры программ

[править]

Консольный рендерер Мандельброта с методом возмущений и предела 10−308 (на языке программирования C++)

[править]

Ключевые особенности

  • Расчёт опорной траектории на 1000 бит всего один раз.
  • Реактивный расчёт миллиарда пикселей на аппаратном double.
  • Расчет можно выполнить существенно быстрее, если использовать билинейную аппроксимацию.
  • При использовании чисел с плавающей запятой двойной точности (порядка 10⁻¹⁵) теория возмущений позволяет приблизиться к уровню 10⁻³⁰⁸ - не дальше.
  • Революционный алгоритм Reference Reset to Zero.
  • Настоящий SSAA 8x8 для идеально сглаженного изображения без алиасинга.
  • Параллелизм OpenMP для высокоскоростного многопоточного рендеринга.
#include <atomic>
#include <cmath>
#include <cstdint>
#include <cstdio>
#include <vector>
#include <algorithm>
#include <mpfr.h>
#include <omp.h>
using std::vector;
using std::min;
const char * CENTER_RE = "-1.99999543561201124623198345433951143502785679245726844745821388800402678499411681518036306219179273434395557574279985918047221291197081186140687781560831995";
const char * CENTER_IM = "-0.00000000000000000000000026198152173811047783694060060607013913873144250985383083459221663448338433592617272786772587281530484110756597337683912309313885172";
const char * VIEW_SIZE = "1.15e-119";
const int WIDTH = 2160;
const int HEIGHT = 2160;
const int AA = 8;
const int MAX_ITER = 50000;
const double ESCAPE_RADIUS_SQUARED = 50000.0;
const int PALETTE_FRAME = 200;
const char * OUTPUT_FILE = "Mandelbrot Set Image 112.bmp";
const mpfr_prec_t PRECISION_BITS = 1000;
const int REF_SIZE = MAX_ITER + 200;
struct Complex {
double re;
double im;
};
#pragma pack(push, 1)
struct BMPHeader {
uint16_t type{0x4D42};
uint32_t size{0};
uint32_t reserved{0};
uint32_t offBits{54};
uint32_t structSize{40};
int32_t width{0};
int32_t height{0};
uint16_t planes{1};
uint16_t bitCount{24};
uint32_t compression{0};
uint32_t sizeImage{0};
int32_t xPixelsPerMeter{2834};
int32_t yPixelsPerMeter{2834};
uint32_t colorsUsed{0};
uint32_t colorsImportant{0};
};
#pragma pack(pop)
int main() {
const double startTime = omp_get_wtime();
const long rawWidth = static_cast<long>(WIDTH) * AA;
const long rawHeight = static_cast<long>(HEIGHT) * AA;
mpfr_t centerRe, centerIm, zReMp, zImMp, tmp1, tmp2, viewSizeMp;
mpfr_inits2(PRECISION_BITS, centerRe, centerIm, zReMp, zImMp, tmp1, tmp2, viewSizeMp, static_cast<mpfr_ptr>(nullptr));
mpfr_set_str(centerRe, CENTER_RE, 10, MPFR_RNDN);
mpfr_set_str(centerIm, CENTER_IM, 10, MPFR_RNDN);
mpfr_set_str(viewSizeMp, VIEW_SIZE, 10, MPFR_RNDN);
const double sampleStep = mpfr_get_d(viewSizeMp, MPFR_RNDN) / rawWidth;
vector<Complex> referenceOrbit;
referenceOrbit.reserve(REF_SIZE);
mpfr_set_ui(zReMp, 0, MPFR_RNDN);
mpfr_set_ui(zImMp, 0, MPFR_RNDN);
for (int iter = 0; iter < REF_SIZE - 1; ++iter) {
Complex z{mpfr_get_d(zReMp, MPFR_RNDN), mpfr_get_d(zImMp, MPFR_RNDN)};
referenceOrbit.push_back(z);
if (z.re * z.re + z.im * z.im > ESCAPE_RADIUS_SQUARED) {
break;
}
mpfr_mul(tmp1, zReMp, zImMp, MPFR_RNDN);
mpfr_sqr(tmp2, zReMp, MPFR_RNDN);
mpfr_sqr(zReMp, zImMp, MPFR_RNDN);
mpfr_sub(zReMp, tmp2, zReMp, MPFR_RNDN);
mpfr_add(zReMp, zReMp, centerRe, MPFR_RNDN);
mpfr_mul_2ui(tmp1, tmp1, 1, MPFR_RNDN);
mpfr_add(zImMp, tmp1, centerIm, MPFR_RNDN);
}
const int referenceLength = static_cast<int>(referenceOrbit.size());
mpfr_clears(centerRe, centerIm, zReMp, zImMp, tmp1, tmp2, viewSizeMp, static_cast<mpfr_ptr>(nullptr));
std::fprintf(stderr, "Reference orbit: %d points\n", referenceLength);
std::fprintf(stderr, "Precomputing skip100 matrices...\n");
vector<Complex> coeff_A(REF_SIZE, {1.0, 0.0});
vector<Complex> coeff_B(REF_SIZE, {0.0, 0.0});
vector<double> rad_R(REF_SIZE, 2.0);
vector<double> aS_squared(REF_SIZE, 0.0);
for (int i = 0; i < referenceLength; ++i) {
double r2 = referenceOrbit[i].re * referenceOrbit[i].re + referenceOrbit[i].im * referenceOrbit[i].im;
aS_squared[i] = (r2 < ESCAPE_RADIUS_SQUARED) ? r2 : 0.0;
}
const int loop_limit = min(static_cast<int>(MAX_ITER), referenceLength - 105);
#pragma omp parallel for
for (int i = 0; i < loop_limit; ++i) {
double min_r2 = ESCAPE_RADIUS_SQUARED;
for (int k = 0; k < 100; ++k) {
if (i + k >= referenceLength) break;
if (aS_squared[i + k] < min_r2) min_r2 = aS_squared[i + k];
}
rad_R[i] = std::sqrt(min_r2);
for (int k = 0; k < 100; ++k) {
if (i + k >= referenceLength) break;
double r_re = referenceOrbit[i + k].re;
double r_im = referenceOrbit[i + k].im;
double next_A_re = 2.0 * (r_re * coeff_A[i].re - r_im * coeff_A[i].im);
double next_A_im = 2.0 * (r_re * coeff_A[i].im + r_im * coeff_A[i].re);
double next_B_re = 2.0 * (r_re * coeff_B[i].re - r_im * coeff_B[i].im) + 1.0;
double next_B_im = 2.0 * (r_re * coeff_B[i].im + r_im * coeff_B[i].re);
coeff_A[i].re = next_A_re; coeff_A[i].im = next_A_im;
coeff_B[i].re = next_B_re; coeff_B[i].im = next_B_im;
}
}
const double PI = 3.14159265358979323846;
uint8_t palette[256][3];
for (int i = 0; i < 255; ++i) {
palette[i][0] = static_cast<uint8_t>(std::lround(127.0 + 127.0 * std::cos(2.0 * PI * i / 255.0)));
palette[i][1] = static_cast<uint8_t>(std::lround(127.0 + 127.0 * std::sin(2.0 * PI * i / 255.0)));
palette[i][2] = palette[i][1];
}
palette[255][0] = 255; palette[255][1] = 255; palette[255][2] = 255;
const int rowBytes = (WIDTH * 3 + 3) & ~3;
vector<uint8_t> image(static_cast<size_t>(rowBytes) * HEIGHT, 0);
std::atomic<int> completedRows{0};
const Complex * reference = referenceOrbit.data();
#pragma omp parallel for schedule(dynamic)
for (int y = 0; y < HEIGHT; ++y) {
uint8_t * row = image.data() + static_cast<size_t>(y) * rowBytes;
for (int x = 0; x < WIDTH; ++x) {
unsigned blueSum = 0; unsigned greenSum = 0; unsigned redSum = 0;
for (int sampleY = 0; sampleY < AA; ++sampleY) {
const double deltaCIm = (static_cast<long>(y) * AA + sampleY - rawHeight / 2) * sampleStep;
for (int sampleX = 0; sampleX < AA; ++sampleX) {
const double deltaCRe = (static_cast<long>(x) * AA + sampleX - rawWidth / 2) * sampleStep;
double deltaRe = 0.0; double deltaIm = 0.0;
double zRe = 0.0; double zIm = 0.0;
int referenceIndex = 0;
int iter = 0;
while (iter < MAX_ITER) {
if (zRe * zRe + zIm * zIm >= ESCAPE_RADIUS_SQUARED) {
break;
}
double eps_abs2 = deltaRe * deltaRe + deltaIm * deltaIm;
double limit_r2 = 1e-60 * rad_R[referenceIndex] * rad_R[referenceIndex];
if (eps_abs2 < limit_r2 && (referenceIndex + 100 < loop_limit) && (iter + 100 < MAX_ITER)) {
double backup_deltaRe = deltaRe; double backup_deltaIm = deltaIm;
int backup_refIdx = referenceIndex; int backup_iter = iter;
double next_eps_re = (coeff_A[referenceIndex].re * deltaRe - coeff_A[referenceIndex].im * deltaIm) + 
(coeff_B[referenceIndex].re * deltaCRe - coeff_B[referenceIndex].im * deltaCIm);
double next_eps_im = (coeff_A[referenceIndex].re * deltaIm + coeff_A[referenceIndex].im * deltaRe) + 
(coeff_B[referenceIndex].re * deltaCIm + coeff_B[referenceIndex].im * deltaCRe);
deltaRe = next_eps_re; deltaIm = next_eps_im;
referenceIndex += 100; iter += 100;
zRe = reference[referenceIndex].re + deltaRe;
zIm = reference[referenceIndex].im + deltaIm;
if (zRe * zRe + zIm * zIm >= ESCAPE_RADIUS_SQUARED) {
deltaRe = backup_deltaRe; deltaIm = backup_deltaIm;
referenceIndex = backup_refIdx; iter = backup_iter;
} else {
continue; 
}
}
const double a = 2.0 * reference[referenceIndex].re + deltaRe;
const double b = 2.0 * reference[referenceIndex].im + deltaIm;
const double nextDeltaRe = a * deltaRe - b * deltaIm + deltaCRe;
deltaIm = a * deltaIm + b * deltaRe + deltaCIm;
deltaRe = nextDeltaRe;
++referenceIndex; ++iter;
zRe = reference[referenceIndex].re + deltaRe;
zIm = reference[referenceIndex].im + deltaIm;
if (zRe * zRe + zIm * zIm < deltaRe * deltaRe + deltaIm * deltaIm || referenceIndex >= loop_limit) {
deltaRe = zRe; deltaIm = zIm; referenceIndex = 0;
}
}
const int remaining = MAX_ITER - iter;
const uint8_t colorIndex = (remaining == 0) ? 255 : static_cast<uint8_t>(remaining % 254);
const int paletteIndex = (colorIndex == 255) ? 255 : (colorIndex - PALETTE_FRAME + 255) % 255;
blueSum += palette[paletteIndex][0];
greenSum += palette[paletteIndex][1];
redSum += palette[paletteIndex][2];
}
}
const int samples = AA * AA;
row[x * 3 + 0] = static_cast<uint8_t>(blueSum / samples);
row[x * 3 + 1] = static_cast<uint8_t>(greenSum / samples);
row[x * 3 + 2] = static_cast<uint8_t>(redSum / samples);
}
const int done = ++completedRows;
if (done % 50 == 0 || done == HEIGHT) {
std::fprintf(stderr, "\rProgress: %d/%d rows (%.1f%%)", done, HEIGHT, 100.0 * done / HEIGHT);
}
}
BMPHeader header;
header.width = WIDTH; header.height = HEIGHT;
header.sizeImage = static_cast<uint32_t>(image.size());
header.size = header.sizeImage + 54;
FILE * file = std::fopen(OUTPUT_FILE, "wb");
if (!file) { std::perror(OUTPUT_FILE); return 1; }
std::fwrite(&header, sizeof(header), 1, file);
std::fwrite(image.data(), 1, image.size(), file);
std::fclose(file);
std::fprintf(stderr, "\nDone! Saved to %s in %.2f seconds.\n", OUTPUT_FILE, omp_get_wtime() - startTime);
return 0;
}


Пример программы построения множества (на языке программирования паскаль)

[править]
uses
  System.Drawing, 
  System.Windows.Forms,
  System.Threading,
  FormsABC; 

procedure DrawMandelbrot(g: Graphics; w,h: integer; scale: real; dx,dy: integer);
const max = 10;
begin
  for var ix:=0 to w-1 do
  for var iy:=0 to h-1 do
  begin
    var x := 0.0;
    var y := 0.0;
    var cx := scale * (ix - dx);
    var cy := scale * (iy - dy);
    var i := 1;
    while i<255 do
    begin
      var x1 := x*x-y*y+cx;
      var y1 := 2*x*y+cy;
      x := x1;
      y := y1;
      if (abs(x)>max) and (abs(y)>max) then break;
      i += 1;
    end;
    if i>=255 then 
      g.FillRectangle(Brushes.Red,ix,iy,1,1)
    else 
      g.FillRectangle(new SolidBrush(Color.FromArgb(255,255-i,255-i)),ix,iy,1,1)
  end;
end;

var 
  Scale := new RealField('Масштаб: ');
  l1 := new FlowBreak;
  dx := new IntegerField('dx: ');
  l2 := new FlowBreak;
  dy := new IntegerField('dy: ');
  l3 := new FlowBreak(20);
  b := new Button(' Нарисовать ');
  p: PaintBox;


procedure Draw;
begin
  var g := p.Graphics;
  DrawMandelbrot(g,p.Width,p.Height,Scale.Value,dx.Value,dy.Value);
  p.Invalidate;
end;

procedure My(o: Object);
begin
  Draw;
end;

procedure Click;
begin
  ThreadPool.QueueUserWorkItem(My);
end;

begin
  MainForm.Title := 'Множество Мандельброта';
  MainForm.SetSize(700, 600);
  MainPanel.Dock := Dockstyle.Left;
  MainPanel.Width := 120;
  Scale.Value := 0.0035;
  dx.Value := 430;
  dy.Value := 280;
  b.Click += Click;

  ParentControl := MainForm;
  p := new PaintBox;
  p.Dock := DockStyle.Fill;
  ThreadPool.QueueUserWorkItem(My);
end.

Пример программы построения множества (на языке программирования C++ с использованием SFML)

[править]
///Heared.hpp
#pragma once

#include <SFML\Graphics.hpp>
//#include <boost\math\bindings\rr.hpp>

typedef long double TYPE;


template<typename T>
struct Complex
{
	T re;
	T im;

	Complex operator += (const Complex & other)
	{
		re += other.re;
		im += other.im;

		return *this;
	}

	Complex operator + (const Complex & other)
	{
		Complex c(*this);
		return c += other;
	}
	Complex operator *= (const Complex & other)
	{
		Complex c;
		c.re = re * other.re - im * other.im;
		c.im = re * other.im + im * other.re;

		*this = c;
		return *this;
	}

	Complex operator * (const Complex & other)
	{
		Complex c(*this);
		return c *= other;
	}
	Complex operator = (const Complex & other)
	{
		re = other.re;
		im = other.im;
		return *this;
	}

	Complex(): re((T)0), im((T)0)
	{
	}

	Complex(const Complex& other) : re(other.re), im(other.im)
	{
	}
	Complex(const T& re, const T& im) : re(re), im(im)
	{
	}
	Complex(const T& re) : re(re), im((T)0)
	{
	}

	T module_sqr()const
	{
		return re * re + im * im;
	}
};



//void drawM(sf::Vector2<TYPE> x, sf::Vector2<TYPE> y);

///draw.cpp
#include "Header.hpp"

void step_clr(sf::Color & clr)
{
	clr.r = clr.r + 5 < 0x100 ? clr.r + 5 : 0xff;
	clr.g = clr.g + 3 < 0x100 ? clr.g + 5 : 0xff;
	clr.b = clr.b + 2 < 0x100 ? clr.b + 5 : 0xff;
}

void drawM(sf::Vector2<TYPE> x, sf::Vector2<TYPE> y)
{
	const TYPE epsilon = 0.005;
	sf::Vector2i sizewindow((int)((x.y - x.x) / epsilon), (int)((y.y - y.x) / epsilon));
	sf::RenderWindow window(sf::VideoMode(sizewindow.x, sizewindow.y), "Mandelbroth");
	sf::Image im;
	const int max_it = 250;
	const int infinity_sqr = 10000;
	im.create(sizewindow.x, sizewindow.y);

	sf::Color clr;
	int xim = 0, yim = 0;
	for (TYPE currx = x.x; currx <= x.y; currx += epsilon, ++xim)
	{
		yim = 0;
		while (xim >= sizewindow.x)
			--xim;
		for (TYPE curry = y.x; curry <= y.y; curry += epsilon, ++yim)
		{
			while (yim >= sizewindow.y)
				--yim;
			clr = sf::Color::Black;
			Complex<TYPE> curr;

			im.setPixel(xim, yim, sf::Color::Black);

			for (int i = 0; i < max_it; ++i)
			{
				if (curr.module_sqr() >= infinity_sqr)
				{
					im.setPixel(xim, yim, clr);
					break;
				}
				curr = curr * curr + Complex<TYPE>(currx, curry);
				step_clr(clr);
			}
		}
	}

	sf::Sprite s;
	sf::Texture t;
	t.loadFromImage(im);
	s.setTexture(t);
	while (window.isOpen())
	{
		sf::Event event;

		while (window.pollEvent(event))

		{

			if (event.type == sf::Event::Closed)

				window.close();

		}

		window.clear();

		window.draw(s);

		window.display();

	}


}

///main.cpp
#include "Header.hpp"
#include "draw.cpp"
#include <iostream>

int main()
{
	TYPE x1, x2, y1, y2;

	//std::cin >> x1 >> x2 >> y1 >> y2;
	x1 = -2.5;
	x2 = 2;

	y1 = -1.5;
	y2 = 1.5;
	drawM(sf::Vector2<TYPE>(x1, x2), sf::Vector2<TYPE>(y1, y2));
}

Добавление цвета

[править]
Фрагмент границы множества Мандельброта в цветном варианте

Строго математически, изображения множеств Мандельброта и Жюлиа должны быть чёрно-белыми. Точка либо попадает внутрь множества, либо нет. Несмотря на это, с помощью компьютера мы можем построить и цветные изображения. Самым распространённым способом является раскрашивание точек снаружи множества в цвет, соответствующий количеству итераций, за которое точка уходит в «бесконечность» или, с точки зрения программы, на определённое расстояние от нуля.

Порядок определения, попадает ли точка z0 внутрь множества (традиционно закрашиваемого чёрным цветом) или нет (закрашивается цветом, зависящим от скорости движения к бесконечности) следующий: на каждой итерации для zn=xn+yn·i вычисляется значение модуля , которое затем сравнивается с «границей бесконечности» (обычно берётся значение, равное 2). Здесь важно обратить внимание, что уже на данном этапе можно ввести определённую оптимизацию вычислений, если проверять не , а , что значительно снизит время расчётов.

Таким образом, если |zn|2 ≤ 4 при любом числе итераций (на практике — при всех вычисленных итерациях), то цвет точки чёрный, в противном случае он зависит от последнего значения n, при котором |zn|2 ≤ 4. Значение n, фактически, обозначает скорость движения zn в бесконечность, и может быть просто индексом в таблице цветов, или использоваться как параметр в более сложном алгоритме.

Данный алгоритм определяет, что если точка удаляется больше чем на 2 от начала координат, то она лежит снаружи множества Мандельброта. Для того, чтобы определить, что точка лежит внутри множества есть много способов. Самое простое решение — ограничить количество итераций неким максимумом. Если точка не вышла за указанную границу, можно считать, что она находится внутри множества.

Точкам около границы множества нужно больше итераций для ухода в бесконечность. Поэтому такие области прорисовываются заметно дольше. Чем дальше от границ множества, тем выше скорость ухода в бесконечность. Для таких точек требуется меньше итераций.

Пример добавления цвета (на PHP)

[править]
<?php

// Default: http://<host>/<dir>/<filename>.php?iter1=64&width=600&height=400&coef=32

function BN($n, $l, $r) {return $n>$l && $n<=$r;}
function SQR($a) {return $a*$a;}

define("COEF",	$_GET["coef"]);
$iter1	=	$_GET["iter1"];
$width	=	$_GET["width"];
$height	=	$_GET["height"];

header("Content-type: image/png");

$img	=	imagecreatetruecolor($width, $height);
$iter2	=	0.01/($width/300);
$yy	=	-1;

for ($y = -1; $y < 1; $y = $y + $iter2) {
	$yy++; $xx=-1;
	for($x = -2; $x < 1; $x = $x + $iter2) {
		$xx++;
		$Cx	=	$x;
		$Cy	=	$y;
		$X	=	$x;
		$Y	=	$y;
		$ix	=	0;
		$iy	=	0;
		$n	=	0;
		while ((SQR($ix) + SQR($iy) < 4) and ($n < $iter1)) {
			$ix 	=	SQR($X) - SQR($Y) + $Cx;
			$iy 	=	2*$X*$Y + $Cy;
			$X	=	$ix;
			$Y	=	$iy;
			$n++;
		}
		if(BN($n,0,7)) $col = imagecolorallocate($img,COEF*$n,0,0);
		elseif(BN($n,7,14)) $col = imagecolorallocate($img,COEF*$n,COEF*$n,0);
		elseif(BN($n,14,21))$col = imagecolorallocate($img,COEF*$n,0,COEF*$n);
		elseif(BN($n,21,28))$col = imagecolorallocate($img,0,COEF*$n,0);
		elseif(BN($n,28,35))$col = imagecolorallocate($img,COEF*$n,COEF*$n,0);
		elseif(BN($n,35,42))$col = imagecolorallocate($img,0,COEF*$n,COEF*$n);
		elseif(BN($n,42,49))$col = imagecolorallocate($img,0,0,COEF*$n);
		elseif(BN($n,49,56))$col = imagecolorallocate($img,COEF*$n,0,COEF*$n);
		elseif(BN($n,56,64))$col = imagecolorallocate($img,0,COEF*$n,COEF*$n);
		imagesetpixel($img, $xx, $yy, $col);
	}
}
imagepng($img);
imagedestroy($img);
?>

Пример добавления цвета (на C++)

[править]
///draw.cpp
#include "Header.hpp"

void step_clr(float & red, float & green, float & blue)
{
	red = (int)(red + 5) < 0xff ? red + 5 : 0xff;
	if((int)red == 0xff)
		green = (int)(green + 3.5) < 0xff ? green + 3.5 : 0xff;
	if((int)green == 0xff)
		blue = (int)(blue + 2.2) < 0xff ? blue + 2.2 : 0xff;
}

void drawM(sf::Vector2<TYPE> x, sf::Vector2<TYPE> y)
{
	const TYPE epsilon = 0.005;
	sf::Vector2i sizewindow((int)((x.y - x.x) / epsilon), (int)((y.y - y.x) / epsilon));
	sf::RenderWindow window(sf::VideoMode(sizewindow.x, sizewindow.y), "Mandelbroth");
	sf::Image im;
	const int max_it = 250;
	const int infinity_sqr = 10000;
	im.create(sizewindow.x, sizewindow.y);

	sf::Color clr;
	int xim = 0, yim = 0;
	for (TYPE currx = x.x; currx <= x.y; currx += epsilon, ++xim)
	{
		yim = 0;
		while (xim >= sizewindow.x)
			--xim;
		for (TYPE curry = y.x; curry <= y.y; curry += epsilon, ++yim)
		{
			while (yim >= sizewindow.y)
				--yim;
			clr = sf::Color::Black;
			Complex<TYPE> curr;
			float red = 0, green = 0, blue = 0;

			im.setPixel(xim, yim, sf::Color::Black);

			for (int i = 0; i < max_it; ++i)
			{
				if (curr.module_sqr() >= infinity_sqr)
				{
					clr = sf::Color(red, green, blue);
					im.setPixel(xim, yim, clr);
					break;
				}
				curr = curr * curr + Complex<TYPE>(currx, curry);
				step_clr(red, green, blue);
			}
		}
	}

	sf::Sprite s;
	sf::Texture t;
	t.loadFromImage(im);
	s.setTexture(t);
	while (window.isOpen())
	{
		sf::Event event;

		while (window.pollEvent(event))

		{

			if (event.type == sf::Event::Closed)

				window.close();

		}

		window.clear();

		window.draw(s);

		window.display();

	}


}
Источник — https://ru.wikibooks.org/w/index.php?title=Множество_Мандельброта&oldid=269747