Я пишу модель Проверку, которая опирается на вычисление коэффициента, который используется интенсивно алгоритмов, заключается в следующем:Как оптимизировать это вычисление
[альтернативный текст] [1]
где q
is double, t
a double and k
a int. e
обозначает экспоненциальную функцию. Этот коэффициент используется на этапах, в которых q
и t
не меняются, а k
всегда начинается с 0 до суммы всех предыдущих коэффициентов (этого шага) достигает 1.
Моя первая реализация была буквальной:
let rec fact k =
match k with
0 | 1 -> 1
| n -> n * (fact (k - 1))
let coeff q t k = exp(-. q *. t) *. ((q *. t) ** (float k)) /. float (fact k)
Конечно, это не так долго, поскольку вычисление всего факториала было просто неосуществимым, когда k
преодолел небольшой порог (15-20): очевидно, результаты начали сходить с ума. Поэтому я перестроил все это делая дополнительных подразделений:
let rec div_by_fact v d =
match d with
1. | 0. -> v
| d -> div_by_fact (v /. d) (d -. 1.)
let coeff q t k = div_by_fact (exp(-. q *. t) *. ((q *. t) ** (float k))) (float k)
Эта версия работает достаточно хорошо, когда q
и t
достаточно «нормальный», но когда что-то становится странно, например, q = 50.0
и t = 100.0
, и я начинаю его вычислить из k = 0 to 100
я получаю серию из 0, за которой следуют NaN с определенного номера до конца.
Конечно, это вызвано операциями с числами, которые начинают слишком близко к 0 или подобным проблемам.
У вас есть идея, как я могу оптимизировать формулу, чтобы иметь возможность давать достаточно точные результаты по широкому распространению входов?
Все должно быть уже 64 бит (поскольку я использую OCaml, который по умолчанию использует удвоение). Возможно, есть способ использовать 128-битные двойники, но я не знаю, как это сделать.
Я использую OCaml, но вы можете предлагать идеи на любом языке, который вы хотите: C, C++, Java и т. Д. Я использовал их все.
Что касается принятия 'log' + http://en.wikipedia.org/wiki/Stirling%27s_approximation – Anycorn
sum (k = 0, ...) x^k/k!== exp (x), поэтому он почти выглядит так, как будто вы делаете exp (-qt) * exp (qt) = 1. Или я что-то упускаю? Заметим, что если qt велико, вы можете использовать exp (qt) = exp (qt/2)^2, т. Е. Вы можете разделить qt на 2 достаточно раза, чтобы сделать серию короткими, а затем квадратизировать это количество раз, чтобы получить желаемый ответ. Не уверен, что это полезно для того, что вы делаете. –