function ds = demo_stiff_withimpl()

    % Vypocetni oblast   
    a = 0.0;
    b = 1.2;
        
    h = 0.006;  pauselen = 1e-8;  ypmax = 2.0;
    %h = 0.02;   pauselen = 0.1;   ypmax = 2.0;
    %h = 0.025;  pauselen = 0.1;   ypmax = 4.0;
    %h = 0.04;   pauselen = 0.1;   ypmax = 6.0;
    %h = 0.12;   pauselen = 0.1;   ypmax = 6.0;
    
    disp(['Krok h nastaven na ' num2str(h)]);
    plot_RK2 = 1;
    
    % Sit ekvidistantnich bodu
    xn=a:h:b;

    % Inicializuj pole pro RK2
    yn2=ones(1,length(xn));
    zn2=zeros(1,length(xn));
    % Zadej pocatecni podminku
    yn2(1)=2;
    zn2(1)=-101;

    % To same pro RK4 a implicitiho Eulera
    yn4=yn2;
    zn4=zn2;
    yn4(1)=yn2(1);
    zn4(1)=zn2(1);

    ynei=yn2;
    znei=zn2;
    ynei(1)=yn2(1);
    znei(1)=zn2(1);

    figure;
    % Nakresli presne reseni
    xp=a:(b-a)/200:b;
    yn=y(xp);
    m=max(max(abs(yn)));
    plot(xp,yn);

    hold on;
    axis([a b 0.0 ypmax]);
   
    % Nakresli numericke reseni: prvni bod (= poc. podm.)
    plot(xn(1),yn2(1),'LineStyle','none','Marker','*','Color',[0 0 0]);
    plot(xn(1),yn4(1),'LineStyle','none','Marker','o','Color',[1 0 0]);
    plot(xn(1),ynei(1),'LineStyle','none','Marker','x','Color',[0 1 0]);	
    legend('funkce','Explic. RK2','Explic. RK4','Impl. Euler','AutoUpdate','off');

    jmenei = 1 + 101*h+100*h^2; 

    for i=1:length(xn)-1,

	pause(pauselen);
        
	% RK2: vypocitej dalsi bod a nakresli ho
	k1 = h*f(xn(i),yn2(i),zn2(i));
	k2 = h*f(xn(i)+h,yn2(i)+k1(1),zn2(i)+k1(2));
	
	yn2(i+1)=yn2(i)+1/2*(k1(1)+k2(1));
	zn2(i+1)=zn2(i)+1/2*(k1(2)+k2(2));
	
        if (plot_RK2==1)
            plot(xn(i+1),yn2(i+1),'LineStyle','none','Marker','*','Color',[0 0 0]);
	end
        
	% RK4: vypocitej dalsi bod a nakresli ho
	k1=h*f(xn(i),yn4(i),zn4(i));
	k2=h*f(xn(i)+1/2*h,yn4(i)+1/2*k1(1),zn4(i)+1/2*k1(2));
	k3=h*f(xn(i)+1/2*h,yn4(i)+1/2*k2(1),zn4(i)+1/2*k2(2));
	k4=h*f(xn(i)+h,yn4(i)+k3(1),zn4(i)+k3(2));
	
	yn4(i+1)=yn4(i)+1/6*(k1(1)+2*k2(1)+2*k3(1)+k4(1));
	zn4(i+1)=zn4(i)+1/6*(k1(2)+2*k2(2)+2*k3(2)+k4(2));
	
	plot(xn(i+1),yn4(i+1),'LineStyle','none','Marker','o','Color',[1 0 0]);

	% Implicitni Euler: vypocitej dalsi bod a nakresli ho
        znei(i+1) = (znei(i) - 100*h*ynei(i))/jmenei;
        %ynei(i+1) = ynei(i) + h*(znei(i) - 100*h*ynei(i))/jmenei;  
        ynei(i+1) = ynei(i) + h*znei(i+1);
  
	plot(xn(i+1),ynei(i+1),'LineStyle','none','Marker','x','Color',[0 1 0]);         
    end;

    hold off;

return;

%%%%%%%%%%%%

function r=f(x,y,z);
%   prava strana   
    r(1)=z;
    r(2)=-100*y-101*z;
    return;

%%%%%%%%%%%%

function r=y(x)
%   presne reseni (jen pro vykreslovani)
    r=exp(-100*x)+exp(-x);
    return;

%%%%%%%%%%%%


