Проект: Мандельброт на Numba и итоги

1. Фрактал: перевод с C++ на Python

Ядро из C08 в Numba-исполнении — сравните с оригиналом, различия только синтаксические (словарь P04 в действии):

import os
# os.environ["NUMBA_ENABLE_CUDASIM"] = "1"   # раскомментируйте на машине без NVIDIA
import numpy as np
from numba import cuda

@cuda.jit
def mandel_kernel(iters, max_iter):
    px, py = cuda.grid(2)                  # двумерный грид — по треду на пиксель
    h, w = iters.shape
    if px >= w or py >= h:
        return                             # защита по обеим осям (C04)

    cr = -2.2 + 3.0 * px / w               # пиксель → точка плоскости
    ci = -1.0 + 2.0 * py / h

    zr = 0.0; zi = 0.0; k = 0
    while k < max_iter and zr*zr + zi*zi <= 4.0:
        zr, zi = zr*zr - zi*zi + cr, 2*zr*zi + ci
        k += 1
    iters[py, px] = k

W, H, MAX_ITER = 900, 600, 200
d_img = cuda.device_array((H, W), dtype=np.int32)       # сразу на устройстве (P05)

block = (16, 16)
grid = ((W + 15) // 16, (H + 15) // 16)
mandel_kernel[grid, block](d_img, MAX_ITER)
img = d_img.copy_to_host()

# сохранение — на этот раз без самодельного PPM: у Python есть PIL
from PIL import Image
t = (img / MAX_ITER) ** 0.45                            # гамма для контраста
rgb = np.dstack([(40 + 200*t), (25 + 170*t), (70 + 120*t)]).astype(np.uint8)
Image.fromarray(rgb).save("mandelbrot.png")

Корректность проверена жёстко: ядро прогнано в симуляторе на уменьшенном поле и сравнено с эталонной последовательной реализацией — совпадение до последнего пикселя (np.array_equal, без допусков: целые числа итераций). Картинка из этого алгоритма — в C08, раздел 3: там она вычислена той же математикой. Замечание о производительности справедливо и здесь, как в C08: у границы множества соседние пиксели требуют разного числа итераций — варпы расходятся, и это видно даже в JS-тренажёре ниже по времени кадра.

2. Случайные числа для тысяч тредов

Монте-Карло из P01 на видеокарте упирается в неожиданное препятствие: обычный генератор случайных чисел — это состояние, которое каждое обращение читает и перезаписывает. Дать тысячам тредов один генератор — устроить гонку (P02) и коррелированные «случайности». Решение CUDA: каждому треду — собственный поток генератора xoroshiro128+, инициализированный так, чтобы потоки не пересекались:

from numba.cuda.random import (create_xoroshiro128p_states,
                                xoroshiro128p_uniform_float32)

@cuda.jit
def mc_pi(states, per_thread, hits):
    t = cuda.grid(1)
    inside = 0
    for i in range(per_thread):
        x = xoroshiro128p_uniform_float32(states, t)   # t — номер СВОЕГО потока
        y = xoroshiro128p_uniform_float32(states, t)
        if x*x + y*y <= 1.0:
            inside += 1
    hits[t] = inside                                   # каждый пишет в свою ячейку — map

threads, blocks, per = 32, 2, 500
states = create_xoroshiro128p_states(threads * blocks, seed=42)
hits = np.zeros(threads * blocks, dtype=np.int32)
mc_pi[blocks, threads](states, per, hits)
print(4.0 * hits.sum() / (threads * blocks * per))    # 3.1172 — проверено

Внутри знакомая механика: каждый тред копит попадания в свой счётчик (map по тредам), финальная сумма — редукция на хосте. Обратите внимание на seed: генератор детерминирован — с одним зерном результат воспроизводится точно, что бесценно для отладки стохастических программ. В симуляторе пример выдаёт ворох предупреждений NumbaWarning при инициализации состояний — это издержки эмуляции, на настоящем GPU их нет.

И обещанный сюрприз из P03: третий target декоратора @vectorize. Формула, написанная для одного элемента, переносится на видеокарту буквально одним словом — Numba сама строит ядро, грид и все копирования данных:

@vectorize(['float32(float32)'], target='cuda')
def expit_gpu(x):
    return 1.0 / (1.0 + math.exp(-x))     # та же сигмоида из P03 — теперь на GPU

3. Итоги траектории: когда C++, когда Numba

Обе линии сошлись в одной точке — один фрактал, две реализации. Пора выдать сводную таблицу выбора:

Практическое правило на выпуск: думайте на языке линии C, прототипируйте на языке линии P. Модель «грид → блоки → варпы», этажи памяти, паттерны map/reduce/scan и критерии рентабельности копирований — общие; какой синтаксис их записывает, вопрос второй. Прототип на Numba за вечер покажет, окупается ли GPU для вашей задачи, — а уже потом решайте, нужна ли пересборка на C++ ради последних процентов и интеграции.

4. Тренажёр: Мандельброт вживую

Контрольные вопросы

• Как проверялась корректность Numba-ядра Мандельброта и почему сравнение без допусков?
  Ядро запускалось в симуляторе на уменьшенном поле и сравнивалось с эталонной последовательной реализацией через np.array_equal. Допусков не нужно: результат — целое число итераций до вылета, а не float; любое расхождение означало бы ошибку логики, а не округления.
• Почему тысячам тредов нельзя пользоваться одним генератором случайных чисел?
  Генератор — это состояние, которое каждое обращение читает и перезаписывает: общий генератор для тредов — гонка данных (испорченное состояние, коррелированные числа). Решение — по независимому потоку генератора на тред: create_xoroshiro128p_states создаёт массив состояний, тред обращается к своему по номеру.
• Где в Монте-Карло-оценке π паттерны map и reduce?
  Map: каждый тред независимо бросает свои per_thread точек и пишет число попаданий в свою ячейку hits[t]. Reduce: суммирование hits — на хосте (их немного) или атомарным сложением/деревом, как в P05. Классическая пара «независимая работа + свёртка результата».
• Что делает @vectorize(target='cuda') и когда он уместнее @cuda.jit?
  Берёт скалярную формулу и сама строит из неё GPU-ядро с гридом и копированиями данных — ни индексов, ни конфигурации запуска. Уместен для чистых поэлементных операций (map). Когда нужны соседи, разделяемая память, атомики или нетривиальная индексация — пишется полноценное @cuda.jit-ядро.
• Сформулируйте правило выбора между линиями C и P для собственного проекта.
  Думать на языке линии C (модель железа, память, паттерны — они общие), прототипировать на линии P: Numba за вечер покажет, окупается ли GPU. Пересборка на C++ оправдана интеграцией в C++-проект, нуждой в специальных API (потоки, cuBLAS, текстуры) или борьбой за последние проценты производительности.

Источники

1. Random Number Generation // Numba Documentation : [сайт]. — URL: https://numba.readthedocs.io/en/latest/cuda/random.html (дата обращения: 09.07.2026).
2. CUDA Ufuncs and Generalized Ufuncs // Numba Documentation : [сайт]. — URL: https://numba.readthedocs.io/en/latest/cuda/ufunc.html (дата обращения: 09.07.2026).
3. CUDA by Numba Examples // NVIDIA Developer Blog : [сайт]. — URL: https://developer.nvidia.com/blog/cuda-by-numba-examples-1/ (дата обращения: 09.07.2026).
4. Mandelbrot set // Wikipedia : [сайт]. — URL: https://en.wikipedia.org/wiki/Mandelbrot_set (дата обращения: 09.07.2026).