Более быстрое сравнение векторов в Юлии - PullRequest
2 голосов
/ 25 июня 2019

Я пытаюсь построить и сравнить, как можно быстрее, два случайных вектора 01 одинаковой длины, используя Юлию, каждый вектор с одинаковым количеством нулей и единиц.

Это все для моделирования Монте-Карло следующего вероятностного вопроса

У нас есть две независимые урны, каждая с n белыми шарами и n черными шарами. Затем мы берем пару шаров, по одному на каждую урну, каждый раз, чтобы опустошить урны. Какова вероятность того, что каждая пара имеет один и тот же цвет?

Я сделал следующее:

using Random
# Auxiliar function that compare the parity, element by element, of two 
# random vectors of length 2n
function comp(n::Int64)
  sum((shuffle!(Vector(1:2*n)) .+ shuffle!(Vector(1:2*n))).%2)
end

Выше генерируются две случайные перестановки вектора от 1 до 2n, добавляется элемент за элементом, применяется модуль 2 к каждому элементу и после суммируются все значения оставшегося вектора. Затем я использую выше четность каждого числа, чтобы смоделировать его цвет: нечетное чёрное и белое чётное.

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

Затем я настраиваю следующую функцию, это просто симуляция МонтеКарло с желаемой вероятностью:

  # Here m is an optional argument that control the amount of random
  # experiments in the simulation
  function sim(n::Int64,m::Int64=24)
  # A counter for the valid cases
  x = 0
  for i in 1:2^m
    # A random pair of vectors is a valid case if they have the
    # the same parity element by element so
   if comp(n) == 0
    x += 1
   end
  end
  # The estimated value
   x/2^m
  end

Теперь я хочу знать, есть ли более быстрый способ сравнения таких векторов. Я попробовал следующую альтернативную конструкцию и сравнение для случайных векторов

shuffle!( repeat([0,1],n)) == shuffle!( repeat([0,1],n))

Затем я изменил соответствующий код на

comp(n)

С этими изменениями код работает немного медленнее, что я тестировал с помощью функции @time. Другими изменениями, которые я сделал, было изменение оператора for для оператора while, но время вычисления осталось прежним.

Поскольку я не программист (на самом деле только вчера я изучил что-то из языка Julia и установил интерфейс Juno), то, вероятно, будет более быстрый способ сделать те же вычисления. Некоторые советы будут оценены, потому что эффективность моделирования MonteCarlo зависит от количества случайных экспериментов, поэтому чем быстрее вычисление, тем большие значения мы можем проверить.

Ответы [ 2 ]

3 голосов
/ 25 июня 2019

Ключевая стоимость в этой задаче составляет shuffle!, поэтому, чтобы максимизировать скорость моделирования, которую вы можете использовать (я добавляю ее в качестве ответа, поскольку она слишком длинна для комментария):

function test(n,m)
  ref = [isodd(i) for i in 1:2n]
  sum(all(view(shuffle!(ref), 1:n)) for i in 1:m) / m
end

Каковы отличия от кода, предложенного в другом ответе:

  1. Вам не нужно shuffle! оба вектора; достаточно shuffle! одного из них, так как результат сравнения инвариантен к любой одинаковой перестановке обоих векторов после их независимого перемешивания; поэтому мы можем предположить, что один вектор после случайной перестановки переставлен, чтобы упорядочить его так, чтобы он имел trues в первых n записях и falses в последних n записях
  2. Я делаю shuffle! на месте (т.е. ref вектор выделяется только один раз)
  3. Я использую функцию all на первой половине вектора; таким образом, проверка останавливается, когда я нажимаю первый false; если я нажму все true в первых n записях, мне не нужно проверять последние n записи, так как я знаю, что они все false, поэтому мне не нужно проверять их
2 голосов
/ 25 июня 2019

Чтобы получить что-то более чистое, вы могли бы генерировать непосредственно векторы значений 0/1, а затем просто позволить Джулии проверять равенство векторов, например,

function rndvec(n::Int64)
 shuffle!(vcat(zeros(Bool,n),ones(Bool,n)))
end

function sim0(n::Int64, m::Int64=24)
  sum(rndvec(n) == rndvec(n) for i in 1:2^m) / 2^m
end

Избегание выделения ускоряет код, как объяснил БогумилКаминьски (и пусть Юлия делает сравнение быстрее, чем его код).

function sim1(n::Int64, m::Int64=24)
  vref = vcat(zeros(Bool,n),ones(Bool,n))
  vshuffled = vref[:]
  sum(shuffle!(vshuffled) == vref for i in 1:2^m) / 2^m
end

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

Я считаю, что это немного не в духе вопроса, но вы могли бы также сделать немного больше математики.Сгенерировано binomial(2*n, n) возможных векторов, и поэтому вы можете просто вычислить

function sim2(n::Int64, m::Int64=24)
   nvec = binomial(2*n, n)
   sum(rand(1:nvec) == 1 for i in 1:2^m) / 2^m
end

Вот некоторые моменты, которые я получаю:

@time show(("sim0", sim0(6, 21)))
@time show(("sim1", sim1(6, 21)))
@time show(("sim2", sim2(6, 21)))
@time test(("test", test(6, 2^21)))

("sim0", 0.0010724067687988281)  4.112159 seconds (12.68 M allocations: 1.131 GiB, 11.47% gc time)
("sim1", 0.0010781288146972656)  0.916075 seconds (19.87 k allocations: 1.092 MiB)
("sim2", 0.0010628700256347656)  0.249432 seconds (23.12 k allocations: 1.258 MiB)
("test", 0.0010166168212890625)  1.180781 seconds (2.14 M allocations: 98.634 MiB, 2.22% gc time)
...