L = 100;
kp = 1;
k = eps;
d = (kp-k)/L;
s = tf('s');
KUS = k;
KOS = 1/(s^3+2*s^2+4*s+2);
K1 = zeros(L,1);
Q = zeros(L,1);
for i=1:1:L
    [y,t] = step(KOS*KUS/(1+KOS*KUS));
    e = 1 - y;
    srednia = mean(e(end-1:end));
    eb = t.*((e - srednia).^2);
    Q(i)= (trapz(t,eb));
    K1(i)=k;
    k = k+d;
    KUS = k;
end
plot(K1,Q);
[min_value,min_index]=min(Q);
disp(K1(min_index));