function X = sylv_tri(T,U,B)
%SYLV_TRI    Solves triangular Sylvester equation.
%            x = SYLV_TRI(T,U,B) solves the Sylvester equation
%            T*X + X*U = B, where T and U are square upper triangular matrices.

m = length(T);
n = length(U);
X = zeros(m,n);

% Forward substitution.
for i = 1:n
    X(:,i) = (T + U(i,i)*eye(m)) \ (B(:,i) - X(:,1:i-1)*U(1:i-1,i));
end