function y = forward(L, b)
%FORWARD. Forward elimination.
% For lower triangular L, y = forward(L, b) solves L*y = b.
    [~, n] = size(L);
    y = zeros(n, 1);
	for k = 1:n
		j = 1:k-1;
		y(k) = (b(k) - L(k,j)*y(j))/L(k,k);
	end
end
