Please check the Code and let me know the error. When i try to run it in MATLAB. It gives me an error .
2 views (last 30 days)
Show older comments
Farooq Zia
on 29 Jun 2017
Commented: Walter Roberson
on 30 Jun 2017
%annual probability that the water height exceeds the threshold hcr.
cv=1
nsam=1000
while cv>0.1
Q=lognrnd(5.1434,0.1416,nsam,1);
Depth=0.2*Q.^0.4329;
Depth1=Depth>4.7;
if any(Depth1)
Pf=sum(Depth1)/nsam
cv=sqrt(1/Pf/nsam);
nsam=round(1/Pf/(0.1^2));
end
end
[SL: indented the while and if blocks to emphasize them.]
1 Comment
Walter Roberson
on 29 Jun 2017
How long should we let it run? It has been more than a minute of computation now but no error message or output has been produced. What error message are you seeing?
Accepted Answer
Walter Roberson
on 29 Jun 2017
It is low probability that even a single element of Q will meet the 0.2*Q.^0.4329 > 4.7 test.
But suppose exactly one element does. Then with nsam = 1000 at first, we have Pf=sum(Depth1)/nsam which would be Pf = 1/1000 . Then cv=sqrt(1/Pf/nsam); would be cv = sqrt(1/(1/1000)/1000) which would be cv = 1 (no change from what it was.) Then nsam=round(1/Pf/(0.1^2)) would be nsam = round(1/(1/1000)/(1/10^2)) which would give nsam = 100000 .
Then we spin for a long time waiting for another item to match. When one finally does, we would have Pf=sum(Depth)/nsam which would be Pf=1/10000 . Then cv=sqrt(1/(1/10000)/10000) would be cv = 1 (no change from what it was). THen nsam=round(1/(1/10000)/(0.1^2) which would take nsam to 10000000 . And so on, continually increasing nsam but getting no closer to reducing cv.
Finally after rather some time, having matched single items N times, nsam = 1000*10^(2*N), we get two matches in depth. Pf=sum(Depth1)/nsam which would be Pf = 2/10^(3+2*N), cv = sqrt(1/Pf/nsam) would be 1/(2/10^(3+2*N)/10^(3+2*N) which is going to give 1/2, and nsam would increase by another factor of 100. But does that help? No, because the next time that a single Depth match is made, cv will be reset back to 1 again.
So in order for cv to reach 1/10, there would have to be 10 or more simultaneous matches on depth out of nsam. What is the probability of that? True, it does increase as nsam increases, but what is our initial probability ?
With lognrnd(5.1434,0.1416,nsam,1)
[M,V] = lognstat(5.1434,0.1416)
M =
173.023129137865
V =
606.311812331174
sqrt(V)
ans =
24.6233996907652
and our target is exp(log(4.7/0.2)/0.4329) == 1469.50114573441 so that is (1469.50114573441 - 173.023129137865) / 24.6233996907652 = 52.6522751885792 standard deviations. That is approximately 8E-605 . That is extremely unlikely to happen in your lifetime.
0 Comments
More Answers (1)
John BG
on 30 Jun 2017
Hi Farooq
1.
increasing Q to 5000
2.
why randomly sampling following a log normal distribution if can enfilade with a straight line?
3.
no need for if any(Depth1) because Depth1 already contains all the values you want to process
cv=1
nsam=1000
while cv>0.1
Q=[1:1:5e3]
Depth=0.2*Q.^0.4329;
Depth1=Depth>4.7;
Pf=sum(Depth1)/nsam
cv=sqrt(1/Pf/nsam);
nsam=round(1/Pf/(0.1^2));
end
Is this of any help?
Pf
= 3.5310
if you find this answer useful would you please be so kind to consider marking my answer as Accepted Answer?
To any other reader, if you find this answer useful please consider clicking on the thumbs-up vote link
thanks in advance
John BG
1 Comment
Walter Roberson
on 30 Jun 2017
John, that will not work. Suppose there is no match the first time, so sum(Depth1) is 0. Then Pf is 0, and 1/Pf is infinity. cv would become infinity. nsam would become infinity as well.
The if any() is important there is important to ensure that there are no divisions by 0 and that cv and nsam are only changed when there is at least one match on depth.
See Also
Categories
Find more on Multidimensional Arrays in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!