forum.vnstele.com
http://forum.vnstele.com/

Обговорення розрахунку формули Пуассона
http://forum.vnstele.com/viewtopic.php?f=9&t=180
Сторінка 1 з 1

Автор:  Mykola Sheremeta [ 23 березня 2012, 00:57 ]
Тема повідомлення:  Обговорення розрахунку формули Пуассона

Мабуть, в десятків студентів під час обчислення імовірності поступлення викликів за допомогою формули Пуассона виникали такі ситуації коли в результаті вони отримували нулі або результати типу NaN(not-a-number) або ще якісь некорректні результати. Це зумовлено недосконалістю використовуваного математичного апарату, а сама недосконалість полягає в тому, що більшість середовищ (MatCAD, MATLAB, Mathematica) коли в результаті виконання якоїсь функції (наприклад функції факторіалу і т.д.) отримують число більше якогось конкретного числа замініють його на значення inf(безмежність, доречі inf це уже символьна змінна), що в подальших обчисленнях приводить до незовсім бажаного результату. (подібне відбувається і з дуже маленькими числами які заміннюються нулем).
Що ж робити в таких ситуаціях?
Слід зробити так щоб программа рахувала не конкретну функцію(наприклад функції факторіалу і т.д.), а весь вираз загалом. Для цього досить часто використовують рекурентну формулу - формулу, що зводить обчислення n- го члена якої-небудь послідовності до обчислення декількох попередніх її членів.
Нижче, я наведу приклад функції яка обчислює формулу Пуассона за допомогою саме такого математичного апарату( функція ця написана в середовищі MATLAB):
Код:
function mass=Poisson(Y,k)
for K=0:k
    for kk=0:K
    if kk==0
        Pk=((Y^kk)*exp(-Y))/ factorial(kk);
    else
        Pk=(Pk*Y)/kk;
    end
    end
    mass(K+1)=Pk;
end
return

Ця функція повертає масив значень імовірності поступлення k викликів, причому елемент масиву з індексом k+1, відповідає імовірності посутплення k викликів (це зумовлено тим, що індексація елментів масиву в середовищі MATLAB починається не з нуля(як в С++), а з одиниці).

Як користуватись цією функцією? Дивіться тут!

Також, я запропоную алогритм для побудови графіка в середовищі MATLAB за допомогою ціє функції:
Код:
at=250/3600; % середня тривалість заняття викликом КС(годин);
N=717; %кількість джерел навантаження;
ac=4.1; %середня кількість викликів від одного джерела за одиницю часу;
k=N; K=0:k;
Y=at*ac*N;
mass=Poisson(Y,k);
graf=plot(K, mass);
set(graf,'lineWidth',2,'Color',[1 0 0]);

title('Графік залежності P(k)(розподіл Пуассона)');
xlabel('k, кількість викликів');
ylabel('P(k), імовірність поступлення k викликів');
grid on;

От власне і все. Дякую, якщо дочитали до кінця. :)

П.С. Якщо ви знайшли в цій функції якусь помилку, або можете запропонувати кращий варіант реалізації цієї функції, пропунуйте свої варіанти, оскільки ця тема призначена саме для цього.
П.П.С. Сподіваюсь викладений тут мною матеріал допоможе вам з виконанням курсової роботи. :D

Автор:  echern [ 23 березня 2012, 15:13 ]
Тема повідомлення:  Re: Обговорення розрахуноку формули Пуассона

Вирішення аналогічної задачі аналогічним способом :D
Код:
% Функція для обчислення імовірності поступлення к викликів (розподіл
% Пуассона. Індекс масиву імовірностей зсунуто на 1, тому що матлаб не
% підтримує індексацію елементів масиву з нуля. (перший елемент масиву -
% для к=0 викликів.
% k - результат виконання функції, дані по осі абсцис (для побудови графіка)
% p - імовірність pk(Y), обчислена за формулою Пуассона
% lambda - інтенсивність поступлення викликів від групи абонентів за одиницю часу.
% t - середня тривалість виклику.
% max_k - Останнє значення К (по осі абсцис)

function [k p] = poissrec (lambda, t, max_k)
Y=lambda*t;
p(1,1)=exp(-Y);
k(1,1)=0;

if (max_k ~= 0)
    for c=1:max_k
        k(c+1, 1) = c;
        p(c+1, 1) = p(c,1)*Y/c;
    end
end
return;


Ну і виклик та перевірка правильності роботи функції:
Код:
% Тестування роботи функції обчислень розподілу Пуассона
% Синьою суцільною лінією - обчислення через рекурентну формулу, червоною
% зірочкою - обчислення безпосередньо за формулою Пуассона.

lambda=10;
t=2;
max_k=50;

[k p] = poissrec (lambda, t, max_k);

for c = 0 : max_k
    pt(c+1, 1) = exp(-lambda*t)*(lambda*t)^c/factorial(c);
end

f=figure;
plot (k, p);
hold on;
plot (k, pt, '*r');

Підписи по осях і сітку на графіку включаємо аналогічно, як в попередньому пості ...
Єдине обмеження по величині Y, це коли e^(-Y) перетворюється в нуль.
Приклад роботи - на рис.1.

Приєднані файли:
test.png
test.png [ 4.77 Кб | Переглянуто 4096 разів ]

Сторінка 1 з 1 Часовий пояс UTC + 2 годин [ DST ]
Powered by phpBB® Forum Software © phpBB Group
https://www.phpbb.com/