function Xs = BisecAllRoots_Schultz(Fun,a,b,TolMax)
solutions = [];
for i = 1:100
    h = (b-a)/(10*i);
    x = a:h:b;
    signs = sign(Fun(x));
    as = [];
    bs = [];
    for i = 1:length(x)-1
        if signs(i+1) ~= signs(i)
            bs(end+1) = x(i+1);
            as(end+1) = x(i);
        end
    end
    if(length(as)~=length(solutions))
        for i = 1:length(as)
            if(length(solutions)>0)
                should_add = 0;
                for j = 1:length(solutions)
                    if solutions(j)>as(i) && solutions(j)<bs(i)
                        %Solution previously found
                    else
                        should_add = 1;
                    end
                end
                if(should_add)
                    solutions(end+1) = BisectionRoot_Schultz(Fun,as(i),bs(i),TolMax);
                end
            else
                solutions(end+1) = BisectionRoot_Schultz(Fun,as(i),bs(i),TolMax);
            end
        end
    else
        disp('No additional Solutions Found')
        Xs=solutions;
        return;
    end
    Xs = 0;
end
end