How do I use the random function in Fsolve Matlab?

1 view (last 30 days)
Hello everyone,
In the following codes, I am trying to determine the variable x that lets beta equal target beta. This does not work for me when I use fsolve. Could anyone please take a look at it and help me with that?
Thank you so much!
% example
clear
clc
% define values
n = 1e4;
target_beta = 2;
x = 2; % predefine x
% parameters
Mean_R = 2000 ; %mean
CoV_R = 0.1;%coefficient of variation
Std_R = Mean_R * CoV_R;
% Monte Carlo Simulation
Ri = Mean_R + Std_R .* norminv(rand(n, 1));
Mean_Q = x * 1100; % Determine which variable, x, will let the beta = target beta
CoV_Q = 0.18;
Std_Q = Mean_Q * CoV_Q;
Qi = Mean_Q + Std_Q .* norminv(rand(n, 1)); % simuate the Q term
g = Ri - Qi; % limit state function g = R - Q
m = find(g < 0); % count the number of failure cases
f = length(m); % find the failure cases number in the g = R-Q vector
pr_failure = f / n; % probability of failure
beta = norminv(1 - pr_failure); %Calculate the reliability index beta
% fsolve
alpha = fsolve(@(x) (beta-target_beta),2)

Accepted Answer

Torsten
Torsten on 22 Mar 2022
This might help to get beta as a function of x:
X = (0*2000/1100:0.1:2000/1100+10).';
n = 1e4;
trials = 1e3;
beta = zeros(numel(X),1);
Pr_failure = zeros(numel(X),1);
for j=1:numel(X)
pr_failure = zeros(trials,1);
x = X(j);
for i = 1:trials
RAND = rand(n,2);
% parameters
Mean_R = 2000 ; %mean
CoV_R = 0.1;%coefficient of variation
Std_R = Mean_R * CoV_R;
% Monte Carlo Simulation
Ri = Mean_R + Std_R .* norminv(RAND(:,1));
Mean_Q = x * 1100; % Determine which variable, x, will let the beta = target beta
CoV_Q = 0.18;
Std_Q = Mean_Q * CoV_Q;
Qi = Mean_Q + Std_Q .* norminv(RAND(:,2)); % simuate the Q term
g = Ri - Qi; % limit state function g = R - Q
m = find(g < 0); % count the number of failure cases
f = length(m); % find the failure cases number in the g = R-Q vector
pr_failure(i) = f / n; % probability of failure
end
pr_failure = mean(pr_failure);
pr_failure = max(eps,pr_failure);
beta(j) = norminv(1 - pr_failure); %Calculate the reliability index beta
Pr_failure(j) = pr_failure;
end
figure(1)
plot(X,beta)
figure(2)
plot(X,Pr_failure)
  4 Comments
Torsten
Torsten on 26 Mar 2022
Are you sure about
beta(j) = norminv(1 - pr_failure);
?
I only know
pr_failure = normcdf(-beta)
thus
beta(j) = - norminv(pr_failure)
Torsten
Torsten on 26 Mar 2022
I was hoping fsolve could be used if the number of trials was big enough to stabilize the distribution for pr_failure (and thus could return stable values to fsolve), but I was not successful. The precision needed to calculate the derivatives in fsolve is too high.

Sign in to comment.

More Answers (1)

Priyanka Kondapalli
Priyanka Kondapalli on 22 Mar 2022
Edited: Priyanka Kondapalli on 24 Mar 2022
Hi,
I do not see any issue with the code provided by you.However, recheck the equation. Please refer to the link below which provides more details on how to use Fsolve.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!