
% function to compute powers a^b  mod n

  function out=modpower(a,b,n)
     out=1;
     for k=1:b
        out=mod(out*a,n);
     end 
  return


