%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                            Michael Pokojovy                             % 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function main;
    figure(1);

    set(gcf, 'PaperUnits', 'centimeters');
    xSize = 18; ySize = 10;
    xLeft = (21 - xSize)/2; yTop = (30 - ySize)/2;
    set(gcf,'PaperPosition', [xLeft yTop xSize ySize]);
    set(gcf,'Position',[0 0 xSize*50 ySize*50]);

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    hold on;
    t = title('Plot der impliziten Funktion $y = y(x)$ mit $e^{y} + y^{3} + 3x^{2} - (e + 4) = 0$, $y(-1) = 1$', 'interpreter', 'latex');
    set(t, 'FontSize', 16);

    F   = @(x, y) exp(y) + y^3 + 3*x^2 - (exp(1) + 4);
    DFy = @(x, y) exp(y) + 3*y^2;

    a = -1 - 3;
    b = -1 + 3;
    
    x = linspace(a, b, 100);
    y = zeros(size(x));
    
    for i = 1:length(x)
        f  = @(y) F(x(i), y);
        df = @(y) DFy(x(i), y);
        
        y(i) = newton(f, df, 1);
    end

    plot(x, y);
    plot(-1, 1, 'ro');
    
    grid on;

    xlabel('x');
    ylabel('y');
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end

function x = newton(f, df, x0)
    max_cnt = 20;
    
    for i = 1:max_cnt
        x = x0 - f(x0)/df(x0);
        
        if (norm(x - x0) < 1E-6) || (norm(f(x)) < 1E-6)
            break;
        else
            x0 = x;
        end
    end
end