function y = mpyH( H, lambda, x )
% mpyH calculate H*x exploiting branch-induced sparsity in H
% mpyH(H,lambda,x) computes H*x where x is a vector and H is a symmetric,
% positive-definite matrix having the property that the nonzero elements on
% row i below the main diagonal appear only in columns lambda(i),
% lambda(lambda(i)), and so on. This is the pattern of branch-induced
% sparsity; and H and lambda can be regarded as the joint-space inertia
% matrix and parent array of a kinematic tree. lambda must satisfy
% 0<=lambda(i)