Метод Монте-Карло (ММ РУиС) — различия между версиями
Материал из SRNS
Korogodin (обсуждение | вклад) (Новая страница: «Навеяно введение в курс Statistical Mechanics Algorithms and Computations with Werner Krauth <source lang="matlab"> % S of square: d^2 % S o…») |
Korogodin (обсуждение | вклад) |
||
Строка 1: | Строка 1: | ||
− | Навеяно | + | Навеяно введением в курс Statistical Mechanics Algorithms and Computations with Werner Krauth |
<source lang="matlab"> | <source lang="matlab"> |
Версия 11:38, 4 февраля 2016
Навеяно введением в курс Statistical Mechanics Algorithms and Computations with Werner Krauth
% S of square: d^2
% S of cicle: pi*r^2 = pi/4 * d^2
% Ss / Sc = 4/pi; pi = Sc / Ss * 4 = Nc / Ns * 4;
N = 10000000;
M = 1000000;
K = N/M;
Nc = 0;
rng('shuffle');
doplot = 0;
for k = 1:K
x = 2*rand(M,1) - 1;
y = 2*rand(M,1) - 1;
inc = find(x.*x + y.*y <= 1);
Nc = Nc + length(inc);
newpi = Nc / (M*k) * 4;
if doplot
plot(x, y, '*b');
hold on
plot(x(inc), y(inc), '*r');
drawnow
end
fprintf('pi = %.6f\n', newpi);
end
% S of cicle: pi*r^2 = pi/4 * d^2
% Ss / Sc = 4/pi; pi = Sc / Ss * 4 = Nc / Ns * 4;
N = 10000000;
M = 1000000;
K = N/M;
Nc = 0;
rng('shuffle');
doplot = 0;
for k = 1:K
x = 2*rand(M,1) - 1;
y = 2*rand(M,1) - 1;
inc = find(x.*x + y.*y <= 1);
Nc = Nc + length(inc);
newpi = Nc / (M*k) * 4;
if doplot
plot(x, y, '*b');
hold on
plot(x(inc), y(inc), '*r');
drawnow
end
fprintf('pi = %.6f\n', newpi);
end