

clc
clear all
close all
warning off

x0 = [  16.6346   11.6442   12.0397   11.2488    5.0218    5.3749   14.5431    6.4924   20.7758   12.4701 ];
lb=  [1;1;1;1;1;1;1;1;1;1]*0.001;
ub=  [20; 20; 20; 20; 20; 20; 20; 20; 20; 20 ]*10;
options = optimset('MaxFunEvals',1e9,'TolFun',1e-15,'TolX',1e-15,'MaxIter',10000,'Display','on');

b=[3.500; 5.000; 2.500; 5.600; 5.600; 5.000; 5.600; 10.500];
h=3.550;
dp=90;
lambda=0.03;
ro=1.12;

b1=b(1);
b2=b(2);
b3=b(3);
b4=b(4);
b5=b(5);
b6=b(6);
b7=b(7);
b8=b(8);

Az=b8*h;
A1=b1*h;
A2=b2*h;
A3=b3*h;
A4=b4*h;
A5=b5*h;
A6=b6*h;
A7=b7*h;
area=[A1 A2 A3 A3 A4 A5 A6 A7 A1 Az];

La=9.580;
Lb=4.700;
Lc=15.030;
Ld=12.070;
Le=8.600;
Lf=8.320;
Lz=5.000;
L=[La Lb Lc Ld Le Lf Lz];

d1=2*b1*h/(b1+h);
d1z=4*b1*h/(2*b1+h);
d2=2*b2*h/(b2+h);
d3=2*b3*h/(b3+h);
d3b=4*b3*h/(2*b3+h);
d4=2*b4*h/(b4+h);
d5=2*b5*h/(b5+h);
d6=2*b6*h/(b6+h);
d6e=2*b6/2*h/(b6/2+h);
d7=2*b7*h/(b7+h);
dz=2*b8*h/(b8+h);
d=[d1 d1z d2 d3 d3b d4 d5 d6 d6e d7 dz];

%ksia=(0.1:1:5);
ksia=0.86;
ksic=0.896;
ksid=1.01;
ksie=0.57;
ksif=0.76;
ksiz=0.14;

for i=1:length(ksia)
    ksia
    ksi=[ksia(i) ksic ksid ksie ksif ksiz];
[x,fval] = lsqnonlin(@(x)potrubi9(x,b,ksi,area,L,dp,lambda,d),x0,lb,ub,options) ;
r(i,:)=x;
hm(i,:)=r(i,:).*area;
hm(i,10)-2*hm(i,8)-hm(i,6)-hm(i,5)
podil(i,:)=hm(i,:)/(hm(i,10)-hm(i,8));
if abs((podil(i,5)+podil(i,6)+podil(i,8))*100-100)>0.1
    x0=x;
    [x,fval] = lsqnonlin(@(x)potrubi9(x,b,ksi,area,L,dp,lambda,d),x0,lb,ub,options) ;
end
x0=x;
end

(podil(:,5)+podil(:,6)+podil(:,8))*100



%figure('name','Prutok na vystupech')
%plot(ksia,podil(:,8)*100)
%hold on
%plot(ksia,podil(:,6)*100,'g')
%plot(ksia,podil(:,5)*100,'r')
%plot(ksia,(podil(:,5)+podil(:,6)+podil(:,8))*100,'m')
%legend('F','D','C','celkem')





%% pro kontrolu rovnic kontinuity%%%%%%
%figure('name','Celkovy prutok')
%plot(ksia,(hm(:,10)-hm(:,8)))

%figure('name','Prutoky usekem A')
%plot(ksia,r(:,1)*A1-r(:,2)*A2)

%figure('name','Prutoky usekem B')
%plot(ksia,r(:,2)*A2-r(:,3)*A3-r(:,4)*A3)

%figure('name','Prutoky usekem C')
%plot(ksia,r(:,4)*A3-r(:,5)*A4)

%figure('name','Prutoky usekem D')
%plot(ksia,r(:,3)*A3-r(:,6)*A5)

%figure('name','Prutoky usekem E')
%plot(ksia,r(:,9)*A1-r(:,7)*A6)

%figure('name','Prutoky usekem F')
%plot(ksia,r(:,7)*A6-2*r(:,8)*A7)

%figure('name','Prutoky usekem Z')
%plot(ksia,r(:,10)*Az-r(:,9)*A1-r(:,1)*A1)


% t=x(10)^2*ro/2*(lambda*Lz/dz)+x(1)^2*ro/2*(ksiz/2)+x(9)^2*ro/2*(ksiz/2);
% w=x(1)^2*ro/2*(lambda*La/2/d1)+x(2)^2*ro/2*(ksia+lambda*La/2/d2);
% e=x(2)^2*ro/2*(lambda*Lb/d3b);%+x(4)^2*ro/2;
% p=x(4)^2*ro/2*(lambda*Lc/2/d3)+x(5)^2*ro/2*(ksic+lambda*Lc/2/d4);
% g=x(3)^2*ro/2*(lambda*Ld/2/d3)+x(6)^2*ro/2*(ksid+lambda*Ld/2/d5);
% z=x(9)^2*ro/2*(lambda*Le/2/d1)+x(7)^2*ro/2*(ksie+lambda*Le/2/d6e);
% l=x(7)^2*ro/2*(lambda*Lf/2/d6)+x(8)^2*ro/2*(ksif+lambda*Lf/2/d7);
% y=[t w e p g z l]


z(1)=(x(10)*Az/2/A1)^2*ro/2*(ksiz);
z(2)=x(2)^2*ro/2*(ksia);
z(3)=x(2)^2*ro/2*(lambda*Lb/d3b);%+x(4)^2*ro/2;
z(4)=x(5)^2*ro/2*(ksic);
z(5)=x(6)^2*ro/2*(ksid);
z(6)=x(7)^2*ro/2*(ksie);
z(7)=x(8)^2*ro/2*(ksif);
zsuma=dp+r(:,10).^2/2*ro-r(:,5).^2/2*ro;
kk=['b','b','b','b','b','b','b']

hh=figure('Position',[50 -200 1600 800])
hold on 
str = {'Z';'A';'B';'C';'D';'E';'F'};
o=[1 2 3 4 5 6 7];
for i=o   
bar(i,z(i),kk(i))
end
%legend('nová geometrie')%
set(gca,'fontsize',20,'fontname','Ariel CE')
title('Úbytek tlaku v jednotlivých úsecích puvodni geometrie')
xlabel('úsek')
ylabel('tlaková ztráta [Pa]')
%ylim([0 100])

set(gca, 'XTickLabel',str, 'XTick',1:numel(str))

grid on
print(hh,'-dtiff','-r300','d:\kluci\kruml\untitled3.tif');



