Генератор случайных чисел с гауссовым распределением для PyOpenCl - PullRequest
Новая
       121

Генератор случайных чисел с гауссовым распределением для PyOpenCl

0 голосов
/ 24 февраля 2020

Генерация гауссовских случайных чисел с использованием numpy оказывается узким местом в моей симуляции Монте-Карло, где я активно использую PyOpenCl. 

np.random.randn(int(1e9))

Поэтому я ищу способ генерировать распределенные по Гауссу случайные числа также с помощью PyopenCl.

Я нашел 6-летнюю тему, задающую аналогичный вопрос. Но я не уверен, как использовать библиотеку VexCL с PyOpenCl: распределенные по Гауссу случайные числа в OpenCL

Любые идеи, как реализовать хороший RNG, который выполняет аналогично np.random.randn в PyOpenCl

1 Ответ

0 голосов
/ 24 февраля 2020

Похоже, что pyopencl уже включает генератор случайных чисел:

https://github.com/inducer/pyopencl/blob/master/pyopencl/clrandom.py

Выполнение простого теста показывает, что среднее и стандартное отклонение сопоставимы с * Реализация 1022 *. Также гистограмма соответствует нормальному распределению с незначительной среднеквадратичной ошибкой.

Кто-нибудь знает дополнительные тесты для проверки качества генератора случайных чисел?

Редактировать: Согласно https://documen.tician.de/pyopencl/array.html 

import pytest
from pyopencl.clrandom import RanluxGenerator, PhiloxGenerator
import pyopencl as cl
import pyopencl.array as cl_array
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm


def test_compare_to_numpy_rand():
    context = cl.create_some_context()
    queue = cl.CommandQueue(context)
    rng =  PhiloxGenerator(context=context)
    mu=0
    sigma=3
    size = (int(1e8))
    dtype = np.float32
    numbers_cl = cl_array.zeros(queue, size, dtype=dtype)
    rng.fill_normal(ary=numbers_cl, mu=mu, sigma=sigma)
    numbers_cl_host = numbers_cl.get()
    mu_meas = np.mean(numbers_cl_host)
    sigma_meas = np.std(numbers_cl_host)

    hist_cl, bins_edges = np.histogram(numbers_cl_host, bins=1000, normed=True)
    delta = bins_edges[1] - bins_edges[0]
    bins_center = bins_edges[:-1]+delta/2
    plt.bar(bins_center, hist_cl,width=delta,label='cl')

    numbers_np = mu + sigma * np.random.randn(size)
    hist_np, bins_edges = np.histogram(numbers_np, bins=bins_edges, normed=True)
    plt.bar(bins_center, hist_np, width=delta,alpha=0.8,label='numpy')

    p = norm.pdf(bins_center, mu, sigma)
    plt.plot(bins_center, p, color='red',label='Exact')
    plt.title('Mean squared error: cl={:.5E} np={:.5E}'.format(np.mean(np.abs(p-hist_cl)**2),np.mean(np.abs(p-hist_np)**2)))
    plt.legend()
    plt.show()

    assert pytest.approx(mu_meas, 1e-2) == mu
    assert pytest.approx(sigma_meas, 1e-2) == sigma
    assert pytest.approx(mu_meas, 1e-2) == numbers_np.mean()
    assert pytest.approx(sigma_meas, 1e-2) == numbers_np.std()

Probability distribution

...