
function [c]=propdamp(m,k,zetas);


[n,n]=size(m);

z=diag(zetas);
[v,d]=eig(inv(m)*k);
mm=v'*m*v;
cn=2*z.*sqrt(d).*mm;
c=inv(v')*cn*inv(v);

