function ps51(a)
if nargin<1, a=1; end
s=tf('s');
sig=hsvd(1/(s+1)+1/(s+a));
p=1/4+1/(4*a);
q=1/(4*a)-1/(a+1)^2;
d=sqrt(p^2-q);
siga=[p+d;p-d];
fprintf('ps51 error for a=%f:  %e  [%f %f]\n', ...
    a,norm(sig-siga),siga(1),siga(2))