clear all


eps=[0 0.00100000000000000	0.00300000000000000	0.00750000000000000	0.0200000000000000	0.0375000000000000	0.0650000000000000	0.100000000000000	0.145000000000000	0.190000000000000	0.218500000000000];
sigma=[0 199665.389473684	349341.620000000	438936.500000000	561585.900000000	653018.700000000		723425.950000000	767033.900000000	791971	799378.800000000	807372.588000000];

fi=0.018;

utau=[0 0.0004 1];
tau=[0 2.87*1000 2.87*1000];


sigmaL=199665;
%epsL=interp1(sigma,eps,sigmaL);
ndelitev=1000;



LzacNabor=[0.80*sqrt(1.5*(sigma(2)/eps(2))*utau(2)*fi/tau(2)):0.005:1.50*sqrt(1.5*(sigma(2)/eps(2))*utau(2)*fi/tau(2))];

%LzacNabor=sqrt(1.5*(sigma(2)/eps(2))*utau(2)*fi/tau(2))
%LzacNabor=[sqrt(1.5*(sigma(2)/eps(2))*utau(2)*fi/tau(2))];


uLnabor=[0.0001:0.00001:0.0005];
tic
for ii=1:size(LzacNabor,2)
    for jj=1:size(uLnabor,2)


Lzac=LzacNabor(ii);
uL=uLnabor(jj);


L=[0:Lzac/ndelitev:Lzac];

uLpoL=L./Lzac*uL;
taupoL=interp1(utau,tau,uLpoL);



sigmapoL=cumsum([0 4/fi*0.5*(taupoL(2:end)+taupoL(1:(end-1))).*diff(L)]);

%sigmapoL(1)=0;
%for i=2:size(L,2)
%sigmapoL(i)=4/fi*((taupoL(i)+taupoL(i-1))/2)*(L(i)-L(i-1))+sigmapoL(i-1);
%end


E=sigma(2)/eps(2);

uLpoLdej=cumsum([0 (0.5*(sigmapoL(2:end)+sigmapoL(1:(end-1)))./E).*diff(L)]);



err(ii,jj)=sigmapoL(end)/sigmaL;
errU(ii,jj)=uLpoLdej(end)/uL;


    end
end
toc

[rr cc1]=find(abs(errU-1)==min(min(abs(errU-1))));
cc=find(abs(err(rr,:)-1)==min(min(abs(err(rr,:)-1))));

%%

% uLpoLdej(1)=0;
% for i=2:size(L,2)
% uLpoLdej(i)=((sigmapoL(i)+sigmapoL(i-1))/(2*E))*(L(i)-L(i-1))+uLpoLdej(i-1);
% end

%Lzac=LzacNabor(1)
%uL=uLnabor(20)

Lzac=LzacNabor(rr);
uL=uLnabor(cc);
L=[0:Lzac/ndelitev:Lzac];

uLpoL=L./Lzac*uL;
taupoL=interp1(utau,tau,uLpoL);

sigmapoL=cumsum([0 4/fi*0.5*(taupoL(2:end)+taupoL(1:(end-1))).*diff(L)]);
E=sigma(2)/eps(2);
uLpoLdej=cumsum([0 (0.5*(sigmapoL(2:end)+sigmapoL(1:(end-1)))./E).*diff(L)]);






figure
hold on
plot(L,uLpoL)
plot(L,uLpoLdej)

figure
hold on
plot(L,taupoL)

figure
hold on
plot(L,sigmapoL)


%%


uLpoL(2)/(L(2)-L(1))

E=sigma(2)/eps(2);




%%

for i=1:size(epspoL,2)-1
   epsPovp(i)=(epspoL(i+1)+epspoL(i))/2;
end
duodX=epsPovp.*diff(L);

uodX(1)=0;
for i=1:size(duodX,2)
 uodX(i+1)=sum(duodX(1:i));   
end







figure
plot(L,uodX)


uL=epsL/2*Lzac;
tauL=interp1(utau,tau,uL);

silaTau=tauL/2*pi*fi*Lzac;
silaSigma=sigmaL*pi*fi^2/4;

silaTau/silaSigma


k=tau(2)/utau(2);
E=sigma(2)/eps(2);

LLL=sqrt(E*fi/k)
