Аппроксимация гамма-функции

Я пытаюсь написать трехпараметрический метод, который аппроксимирует гамма-функцию на определенном интервале. Аппроксимация должна быть суммой Римана с правым концом.

Гамма-функция определяется как:

GAMMA(s) = 
inf
INT  x^(s-1) * exp(-x) dx
0

Таким образом, аппроксимация суммы Римана на правом конце в интервале (0, m) должна быть:

GAMMA(s) ~  
m
SUM  ((m/n)*i)^(s-1) * exp(-(m/n)*i) * delta_x        where delta_x = (m/n)
i=1

Мой код выглядит следующим образом:

def gamma(x = 4.0, n = 100000, m = 2500)
  array = *(1..n)
  result = array.inject(0) {|sum, i| sum + ((((m/n)*i)**(x-1))*((2.7183)**(-(m/n)*i))*(m/n))}  
end

puts gamma

Код должен возвращать приближение для 3! = 6, но вместо этого возвращает 0,0. Есть идеи, где я могу ошибаться?


ruby gamma-function
person Ice101781    schedule 15.08.2014    source источник

Ответы (2)


arrow_upward
3
arrow_downward

Проблема в том, что когда вы делаете m/n

вы делаете целочисленное деление (например, 3/4 = 0), когда вы ожидаете деление с плавающей запятой (3/4 = 0,75)

вам нужно определить ваши n и m как поплавки.

Вы можете переписать его как

def gamma(x = 4.0, n = 100000, m = 2500)
  n = n.to_f
  m = m.to_f

  (1..n).to_a.inject(0) do |sum, i|
    sum + ((((m/n)*i)**(x-1))*((Math::E)**(-(m/n)*i))*(m/n))
  end  
end

PS: также вам не нужны переменные array и result. PS2: рассмотрите возможность использования Math::E вместо 2.7183

person xlembouras    schedule 15.08.2014
comment
Я только что реализовал процесс в Excel и получил правильный ответ, поэтому я знал, что мой подход был правильным. Я кодировщик n00b, поэтому я решил, что это должно быть так. Спасибо еще раз. Этот сайт отличный. - person Ice101781; 16.08.2014

arrow_upward
2
arrow_downward

Ваша проблема была определена @xlembouras. Вы можете написать метод следующим образом.

Код

def gamma(x = 4.0, n = 100000, m = 2500)
  ratio = m.to_f/n
  xm1 = x-1.0
  ratio * (1..m).inject(0) do |sum,i|
    ixratio = i*ratio
    sum + ixratio**xm1 * Math.exp(-ixratio)
  end
end

Примеры

gamma(x=4.0, n= 40, m =10).round(6) #=> 1.616233
gamma.round(6)                      #=> 6.0

Пожалуйста, подтвердите, что эти расчеты верны.

person Cary Swoveland    schedule 15.08.2014
comment
гамма(4) должна равняться 3! = 6 - person Ice101781; 16.08.2014
comment
Я планирую написать метод, который расширит функцию Гамма на всю комплексную плоскость. Конечная цель состоит в том, чтобы затем написать метод, который позволит мне оценивать значения функционального уравнения дзета-функции Римана. - person Ice101781; 16.08.2014