Problem in Monte Carlo simulation.
조회 수: 1 (최근 30일)
이전 댓글 표시
My code is not giving me the desired result. It always shows 1. For your ease let me explain the base of code again
The failure occur when "C >= Cr". so I intend to calculate the probabilty of that failure, lets say its P_fail. Here, C will be Calculated using
C(i,1) = Cs(i,1).*(1-erf(Cd(i,1)./(2.*sqrt(con))));
and Cr = 1.2
and ultimately final output is calculated when P_fail becomes greater than P_fail_Max (which equals to 0.10).
clear all
clc
n = 200;
%cover depth
mean_Cd = 16.1;
sigma_Cd = 0.23;
%diffusion coefficient at 28 days
mean_D = 3.87*10^(-12);
sigma_D = 0.52;
%surface chloride concn
mean_Cs = 13.1;
sigma_Cs = 0.103;
%Time exponent
mean_T = 0.2;
sigma_T = 0.2;
%critical chloride content
mean_Cr = 1.2;
sigma_Cr = 0.2;
%Now use normrnd function and lognrnd to generate values randomly
Cd = normrnd(mean_Cd, sigma_Cd, [n,1]);
D = lognrnd(mean_D, sigma_D, [n,1]);
Cs = normrnd(mean_Cs, sigma_Cs, [n,1]);
T = normrnd(mean_T, sigma_T, [n,1]);
Cr = normrnd(mean_Cr, sigma_Cr, [n,1]);
%Create for loop
for i=1:n
rng(0,'twister'); %generating a vector of 200 random values for Cr values
a = 0.52;
b = 3.87*10^(-12);
D = a.*randn(200,1) + b;
%display(D)
rng(0,'twister'); %generating a vector of 200 random values for Cr values
a = 0.2;
b = 0.2;
T = a.*randn(200,1) + b;
% display(T)
if (D(i,1).*T(i,1) > 0 ) % Check to see if this number is real and positive before calculating sqrt
con = D(i,1).*T(i,1);
else
con = 1;
%break;
end
rng(0,'twister'); %generating a vector of 200 random values for Cr values
a = 0.2;
b = 1.2;
Cr = a.*randn(200,1) + b;
%display(Cr)
rng(0,'twister'); %generating a vector of 200 random values for Cr values
a = 0.23;
b = 16.1;
Cd = a.*randn(200,1) + b;
% display(Cd)
rng(0,'twister'); %generating a vector of 1000 random values for Cr values
a = 0.103;
b = 13.1;
Cs = a.*randn(200,1) + b;
%display(Cs)
% MAIN EQUATION TO CALCULATE C
C(i,1) = Cs(i,1).*(1-erf(Cd(i,1)./(2.*sqrt(con))));
a = 1;
b = i;
C = a.*randn(20,1) + b;
display(C)
Pf = mean(C > Cr);
display(Pf)
end
Pfmax = 0.10;
Tser = mean(Pf >= Pfmax);
display(Tser) % final output
댓글 수: 2
Walter Roberson
2022년 5월 27일
https://www.mathworks.com/matlabcentral/answers/1707460-how-to-code-monte-carlo-simulation#comment_2181870
답변 (0개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Signal Generation and Preprocessing에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!