% 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

%---------------------------------

% function to find inverse mod n
  function out=modinverse(a,n);
     x=mod(a*[1:n],n);
     out=find(x==1);
  return;

%-------------------------------

 % Carmichael's totient function 

   function out=lambda(p,q)
     out=lcm(p-1,q-1)
   return

%------------------------------

%   time_factor.m  measure time needed to factor numbers
N=50;
  t=zeros(1,N);
  for k=1:N
     for j=1:30
        tic;
        factor(randi(2^k));
        t(k)=toc+t(k);
      end
  end
  plot(log(t/N))