function [mc, ec, Lckm, Lkm]=spl_f(n,x,xpl,mtd,mnl,mpl,Ftd,Fnl,Fpl,delta,pp)
%vystupni parametry: Lckm - matice odlehlosti pro mc, ec; Lkm - matice
%odlehlosti pro vesaky
%funkce pro vypocet elasticity a kmitajici hmoty vedeni s pridavakem pres
%jeden vesak

F=100;													%samozrejme parametry
g=9.81;

%vytvoreni hodnot Ftd/x a celk. matice elementarnich vah pro TD
%prvni prvek matic odpovida prvnimu vesaku definovaneho kotevniho useku,
%posledni prvek matice xtdc odpovida predposlednimu vesaku def. kot. useku

c=1;
for a=1:pp
    xtdc(c)=x(a,2);                         %odlehlost za prvnim vesakem
    Mtdc(c)=(x(a,1)+x(a,2)/2)*mtd+(xpl(a,1)-x(a,1))/2*mpl+x(a,1)*mpl;  %prvni vesak drzi vahu troleje a PL uz od bocaku
    if a==1
        Lkm(c)=x(a,1);                      %kumulativni odlehlost uplne prvniho vesaku
    else
        Lkm(c)=Lkm(c-1)+x(a-1,n(a-1)-1)+x(a,1);      %kumulativni odlehlost prvniho vesaku dalsich poli zabira i delku od posledniho vesaku predch. pole k bocaku
    end
    c=c+1;  
    for b=3:n(a)-2                           %pro druhy az predposledni vesak
        xtdc(c)=x(a,b);                      %odlehlost do matice x
        Mtdc(c)=(x(a,b-1)+x(a,b))/2*mtd;     %vaha pro TD
        Lkm(c)=Lkm(c-1)+x(a,b-1);            %kumulativni odlehlost
        c=c+1;
    end
    Mtdc(c)=(x(a,n(a)-2)/2+x(a,n(a)-1))*mtd+(xpl(a,2)-x(a,n(a)-1))/2*mpl+x(a,n(a)-1)*mpl; %posledni vesak drzi vahu troleje az k bocaku
    if a~=pp                                 %odlehlost nelze pro uplne posledni vesak
        xtdc(c)=x(a,n(a)-1)+x(a+1,1);        %odlehlost posl. vesaku zapocitava i vzdalenost od bocaku do prvniho vesaku nasl. rozpeti
        Lkm(c)=Lkm(c-1)+x(a,n(a)-2);         %kumulativni odlehlost klasicky
    end
    c=c+1;
end
ntd=c-2;                                     %pro posledni vesak nemame odlehlost - nepocitame
c=1;

%vytvoreni hodnot Fnl/x a celk. matice elementarnich vah pro NL
%prvni prvek matic odpovida prvnimu vesaku definovaneho kotevniho useku,
%posledni prvek matic odpovida poslednimu vesaku def. kot. useku,
%i zavesy maji svuj prvek: length(znlc)=length(xnlc)=length(ztdc)+pp-1
%vektor fp jednotkovou hodnotou prvku upozornuje na index, ktery odpovida
%prave zavesu

fp=zeros(1,pp*16);
for a=1:pp
    if a~=1                     %pro zaves krome prvniho
        xnlc(c)=xpl(a,1);       %pro zaves pocitame pouze odlehlost, jeho vaha se nehybe
        fp(c)=1;                %indikace zavesu
        c=c+1;
    end
        xnlc(c)=x(a,1)+x(a,2)-xpl(a,1);                 %pripojeni pridavaku
        Mnlc(c)=(x(a,1)+x(a,2))/2*mnl+(xpl(a,1)-x(a,1))/2*mpl;   %pripocitame i kousek PL
        c=c+1;
        xnlc(c)=x(a,3);                                 %prvni vesak na NL
        Mnlc(c)=(x(a,1)+x(a,2)-xpl(a,1)+x(a,3))/2*mnl;  %vaha pro prvni vesak na NL
        c=c+1;
    for b=4:n(a)-3                                      %pro druhy az predposledni vesak
        xnlc(c)=x(a,b);                                 %odlehlost nasl. vesaku
        Mnlc(c)=(x(a,b-1)+x(a,b))/2*mnl;                %vaha odpovidajic vesaku
        c=c+1;
    end
        xnlc(c)=x(a,n(a)-2)+x(a,n(a)-1)-xpl(a,2);       %posledni vesak
        Mnlc(c)=(x(a,n(a)-3)+x(a,n(a)-2)+x(a,n(a)-1)-xpl(a,2))/2*mnl;
        c=c+1;                                        
        xnlc(c)=xpl(a,2);                               %pripojeni pridavaku
        Mnlc(c)=(x(a,n(a)-2)+x(a,n(a)-1))/2*mnl+(xpl(a,2)-x(a,n(a)-1))/2*mpl; %pripocitame i kousek PL
        c=c+1;
end
nnl=c-2;
xnlc(c-1)=[];


%vytvoreni celkove matice pro nosne lano x,vyskrtnuti sloupcu a radku, ktere odpovidaji
%"hmote" v zavesu, pridani dvou prazdnych radku pro udaje prvniho vesaku a pridani
%rovnic pro pridav. lano. Vznikne matice, ze ktere se mohou primo brat podmatice
%14x14 pro vlastni vypocet. Pro trolejovy drat celkovou matici x nelze
%vytvorit, nebot jeji prvky jsou meneny dosazovanim pritlacne sily v miste promenne 
%odlehlosti lambda 

mxnlc=zeros(nnl);                       %na kazdy zaves se prida jeden radek a sloupec  
for a=2:nnl                                   %naplneni matice hodnotama pro nl
    mxnlc(a,a-1)=Fnl/xnlc(a-1);
    mxnlc(a,a)=-Fnl*(1/xnlc(a-1)+1/xnlc(a));
    mxnlc(a,a+1)=Fnl/xnlc(a);
end
a=1;
%pom=zeros(1,length(Mnlc)+1);           %doplneni uvodni nuly do matice Mnlc
%pom(2:end)=Mnlc;
%Mnlc=pom;

fnl=zeros(1,pp*17);
ftd=zeros(1,pp*17);
app=0;
while a<nnl
    if fp(a)==1                     %pokud je zaves, musime doplnit prvky na pridavak!
      %rozsireni matice mxnlc      
      d=length(mxnlc);                          %budeme vkladat radky do matice mxnlc
      mxnlc(:,d+1)=0;                           %matici rozsirime
      mxnlc(d+1,:)=0;
      mxnlc(a+2:d+1,a+2:d+1)=mxnlc(a+1:d,a+1:d);%odsuneme prave dolni prvky 
      mxnlc(:,a:a+1)=0;                         %vynulujeme vznikle pasy
      mxnlc(a:a+1,:)=0;
      %rozsireni vektoru Mnlc
      Mnlc(length(Mnlc)+1)=0;                   %prodlouzime matici vah
      Mnlc(a+2:end)=Mnlc(a+1:end-1);            %odsuneme konec dal
      Mnlc(a+1)=0;                              %vynulujeme prvek
      %vytvoreni pomocne matice pro pridavak. Pokud jsme jiz prosli nejake
      %zavesy (app), jsme v matici mxnlc o app prvku dal nez v xnlc a o
      %2*app prvku dal nez v xtdc
      pom=zeros(4);
      pom(1,1)=-Fpl/(xtdc(a-2-2*app)-xnlc(a-2-app));
      pom(1,2)=-pom(1,1);
      pom(2,1)=pom(1,2);
      pom(2,2)=-Fpl/(xtdc(a-2-2*app)-xnlc(a-2-app))-Fpl/xtdc(a-1-2*app);
      pom(2,3)=Fpl/xtdc(a-1-2*app);
      pom(3,2)=pom(2,3);
      pom(3,3)=-Fpl/xtdc(a-1-2*app)-Fpl/(xtdc(a-2*app)-xnlc(a+1-app));
      pom(3,4)=Fpl/(xtdc(a-2*app)-xnlc(a+1-app));
      pom(4,3)=pom(3,4);
      pom(4,4)=-pom(3,4);
      mxnlc(a-1:a+2,a-1:a+2)=mxnlc(a-1:a+2,a-1:a+2)+pom;
      %zbyle informacni vypocty
      fp(a+2:length(fp))=fp(a+1:length(fp)-1);  %posuneme vektor fp
      fnl(a-1)=1;                               %indikace radku vyhradne pro NL
      fnl(a+2)=1;
      ftd(a)=1;                                 %indikace radku vyhradne pro TD
      ftd(a+1)=1;
      nnl=nnl+1;
      app=app+1;                                %posunuti na dalsi pole
    end
    a=a+1;                                      %na dalsi vesak
end

%v tuto chvili jsou predpripravene prvky do matice mxtd (xtdc) a celkova matice 
%pro mxnl (mxnlc). Je treba pomoci promenne odlehlosti vytvorit matici
%mxtd, matici mxnl rozsirit o jeden radek pro sberac a matici mxtd rozsirit
%o radky a sloupce, ktere odpovidaji zavesum pridavaku. Obe matice lze pak
%secist a vypocitat zdvih.
pom=1;
for a=1:nnl
    if fnl(a)==0                                %vytvoreni vektoru odpovidajicich prvku
        ftd(pom)=a;                             %kazdemu prvku matice mxtd se priradi prvek
        pom=pom+1;                              %z matice mxnl
    end
end

c=1;
for a=1:ntd-11                                  %pro vsechny lanicka
    for L=delta/2:delta:xtdc(a+5)-delta/2       %podle delkoveho kroku "delta"
        Lckm(c)=Lkm(a+5)+L;                     %kumulativni odlehlosti mist vypoctu vypoctena pomoci kumulativnich odlehlosti vesaku. Zaciname od sesteho lanicka.

        %vytvoreni matice vypoctove TD a NL
        
        pomxtd(1:5)=xtdc(a:a+4);
        pomxtd(6)=L;
        pomxtd(7)=xtdc(a+5)-L;
        pomxtd(8:12)=xtdc(a+6:a+10);
        pomMtd(1:5)=Mtdc(a+1:a+5);
        pomMtd(6)=0;
        pomMtd(7:11)=Mtdc(a+6:a+10);
        for b=2:12                              %vytvoreni matice pro TD do vypoctu
            mxtd(b,b-1)=1/pomxtd(b-1);
            mxtd(b,b)=-1/pomxtd(b-1)-1/pomxtd(b);
            mxtd(b,b+1)=1/pomxtd(b);
        end
        mxtd=mxtd.*Ftd;
        mxtd=mxtd(2:12,2:12);                   %orezani matice
        pomftd(1:5)=ftd(a+1:a+5);               %vytvoreni vektoru odpovidajicich prvku do 
        pomftd(7:11)=ftd(a+6:a+10);             %vypoctu
        pomftd(6)=ftd(a+5);
        
        z=zeros(11,1);                          %prava strana
        z(6)=-F;
        
        pomMnl=Mnlc(pomftd(1):pomftd(11));      %matice vah NL do vypoctu
        pommxnl=mxnlc(pomftd(1):pomftd(11),pomftd(1):pomftd(11)); %matice pro NL do vypoctu
        n=length(pommxnl)+1;                     %takova bude celkova velikost vypoctove matice
        
        %rozsirovani matice NL
        pommxnl(:,pomftd(6)-pomftd(1)+3:n)=pommxnl(:,pomftd(6)-pomftd(1)+2:n-1);
        pommxnl(pomftd(6)-pomftd(1)+3:n,:)=pommxnl(pomftd(6)-pomftd(1)+2:n-1,:);
        pommxnl(:,pomftd(6)-pomftd(1)+2)=zeros(n,1);
        pommxnl(pomftd(6)-pomftd(1)+2,:)=zeros(1,n);
        pomMnl(pomftd(6)-pomftd(1)+3:n)=pomMnl(pomftd(6)-pomftd(1)+2:n-1);
        pomMnl(pomftd(6)-pomftd(1)+2)=0;
        
        %rozsirovani matice TD na pozicich uchyceni pridavaku
        for b=2:n
            if pomftd(b)-pomftd(b-1)>1
                mxtd(:,b+1:end+1)=mxtd(:,b:end);    %rozsireni vypoctove matice TD
                mxtd(b+1:end+1,:)=mxtd(b:end,:);
                mxtd(:,b)=zeros(length(mxtd),1);
                mxtd(b,:)=zeros(1,length(mxtd));
                pomMtd(b+1:end+1)=pomMtd(b:end);    %rozsireni matice vah TD
                pomMtd(b)=0;
                pomftd(b+1:end+1)=pomftd(b:end);    %rozsireni vektoru odpovidajic.prvku
                pomftd(b)=pomftd(b)-1;
                z(b+1:end+1)=z(b:end);              %rozsireni vektoru prave strany
                z(b)=0;
            end
        end
        mx=pommxnl+mxtd;
        pomM=pomMtd+pomMnl;
        y=mx\z;
        ec(c)=max(y);
        pomM=pomM*y;
        mc(c)=pomM/ec(c);
        ec(c)=ec(c)/F;
        c=c+1;
        clear pomMtd;
        clear pomMnl;
        clear pommxnl;
        clear mxtd;
        clear pomftd;
    end
end

       