matlab 蒙特卡罗仿真

发布网友发布时间：2022-05-01 09:43

我来回答

共2个回答

热心网友时间：2022-06-27 18:24

之前回答过题主的另一个相关问题（编号1384750375215298220），其中的随机数生成稍微有点问题，请把其中的

r = rand(N,1)*D/2;

改为

r = sqrt(rand)*D/2; % 半径的概率与其长度一致

以使得随机数在圆内均匀分布。

这个问题的改动主要是考虑了守门员这个因素，另外图形显示的实现上也有些变化。参考代码如下：

% 参数输入
d = inputdlg('射门次数R','试验设置',1,{'100'});
if isempty(d), return, end
R = round(str2double(d{1}));
% R = 100;
L = 4;
W = 2;
D = sqrt(L^2+W^2);

% 绘图
clf
t = linspace(0,2*pi,200);
plot(D/2*cos(t),D/2*sin(t),'linewidth',2);
hold on
[x,y] = meshgrid((-2:2)*L/4, (-1:1)*W/2);
surf(x,y,x*0,'Facealpha',0.3,'Edgealpha',0.5)
h = plot(NaN,NaN,'ro');
k1 = patch(NaN,NaN,'k','Facealpha',0.6);
k2 = patch(NaN,NaN,'k','Facealpha',0.6);
axis equal

% 守门员区域定义
K{1} = [-L/4 L/4 -W/2 W/2];
K{2} = [-L/2 -L/4 0 W/2; -L/4 0 -W/2 0];
K{3} = [0 L/4 -W/2 0; L/4 L/2 0 W/2];
K{4} = [-L/2 0 -W/2 0];
K{5} = [0 L/2 -W/2 0];

% 模拟
P = zeros(R,1);
Q = P;
for n = 1 : R
    r = sqrt(rand)*D/2;  % 半径的概率与其长度一致
    t = rand*2*pi;
    x = r.*cos(t);
    y = r.*sin(t);

    % 判断在球门范围内
    P(n) = abs(x) <= L/2 & abs(y) <= W/2;

    % 判断是否被守门员扑出
    k = ceil(rand*5);
    if any( x>K{k}(:,1) & x<K{k}(:,2) & y>K{k}(:,3) & y<K{k}(:,4) )
        Q(n) = 1;
    end

    % 更新绘图
    set(k1,'xData',K{k}(1,[1 2 2 1]),'yData',K{k}(1,[3 3 4 4]))
    if size(K{k},1) > 1
        set(k2,'xData',K{k}(2,[1 2 2 1]),'yData',K{k}(2,[3 3 4 4]))
    else
        set(k2,'xData',NaN,'yData',NaN)
    end

    set(h,'xdata',x,'ydata',y);
    if ~P(n)
        set(h, 'color','m');
    elseif P(n) && ~Q(n)
        set(h, 'color','g');
    else
        set(h, 'color','r');
    end
    drawnow
end

% 计算概率
p = mean(P)
q = mean(Q)
mean(P.*~Q)

热心网友时间：2022-06-27 18:25

贴一个蒙特卡洛方法的matlab程序，供大家使用。

{3 x& K/ i1 i( D8 C0 c$ O
% Example Monte Carlo Simulation in Matlab 0 O5 \; P" t# t7 v8 c& @
% Function: y = x2^2/x1 5 Z0 W4 e9 q, d5 B+ c
%
% Generate n samples from a normal distribution 4 s! c6 y, I6 H" d) K+ v. Y; X: Q
% r = ( randn(n,1) * sd ) + mu 4 U F* Q) t, T# q* w/ K' Q
% mu : mean / E( P8 U" c* o! G8 s/ x
% sd : standard deviation
%
% Generate n samples from a uniform distribution 2 u# ^& K. [0 z% F) @1 y
% r = a + rand(n,1) * (b-a) - D+ }& U$ w- M9 @& Q9 W, Z
% a : minimum
% b : maximum
n = 100000; % The number of function evaluations 7 x5 a" @- F& O- Z; w5 j
% --- Generate vectors of random inputs ! K& x0 ^# X+ q( V6 {
% x1 ~ Normal distribution N(mean=100,sd=5)
% x2 ~ Uniform distribution U(a=5,b=15)
x1 = ( randn(n,1) * 5 ) + 100; 2 B' l3 n) V) D$ ~
x2 = 5 + rand(n,1) * ( 15 - 5 ); \: O: Y( w3 [9 d: V4 r( k4 {
% --- Run the simulation
% Note the use of element-wise multiplication - ~% x$ `7 A6 v9 R* F
y = x2.^2 ./ x1; ' g$ O7 U; R* F% `
% --- Create a histogram of the results (50 bins)
hist(y,50); / M9 m+ s( [* w" J2 I% s/ X
% --- Calculate summary statistics
y_mean = mean(y)
y_std = std(y) ; R7 A2 y M/ T" p, h* m
y_median = median(y)