function u = tridproduct(ld, d, ud, x)
% Multiply the tridiagonal matrix A = [ld d ud] by x
% and return the product in u
% Input:
%     ld - lower diagonal
%     d  - diagonal
%     ud - upper diagonal
%     x  - rhs of matrix product
% Output:
%     u - A*x
 
n = length(d);
u = zeros(n,1);
u(1) = d(1)*x(1) + ud(1)*x(2);
u(n) = ld(n-1)*x(n-1) + d(n)*x(n);
for i=2:n-1
  u(i) = ld(i-1)*x(i-1) + d(i)*x(i) + ud(i)*x(i+1);
end
