KnuthS(v,n)={
  my(u=vector(n,i,i));
  for(i=n+1,#v,
    if(random(i)<n,u[random(n)+1]=i)
  );
  vecextract(v,u)
};
test()={
  my(v=vector(10),t);
  for(i=1,1e5,
    t=KnuthS([0,1,2,3,4,5,6,7,8,9],3);
    v[t[1]+1]++;v[t[2]+1]++;v[t[3]+1]++
  );
  v
};