An example
clc,clear
% umaxa=163;
umaxa=113;
umaxb=167;
n = 20;
area=0.064;
rng('shuffle')
Umax=umaxa + (umaxb-umaxa).*rand(n,1);
nn = n;
cm = jet(100); % colormap
F = 1; % start force
while (nn > 0) && (F < 10000)
stress = (F*50)/(nn*area); % stress of each fiber
on = stress < Umax; % fibers remain
nn = sum(on); % actual number of fibres
% F = F + 1;
F = F + nn/n/10;
cla
title(sprintf('Stress = %0.2f',stress))
hold on
for i = 1:n
if on(i) % if fiber exists
ind = round(100*stress/Umax(i)); % color index
plot([i i],[0 n],'color',cm(ind,:)) % plot fiber
str = sprintf('Umax = %0.0f',Umax(i)); % label
text(i,n/2,str,'ROtation',90) % plot label
end
end
hold off
pause(1)
end