clear
global Kop
etts=[];
multtimes=[];
n1s=[];
enorms=[];
for ii=1:5,
n1=39+2*ii;
tn1=n1*(n1+1)/2;
m=2*n1^2;n=2*tn1;
%rhs=randn(m,1)+.01*ones(m,1);
rhs=sprandn(m,1,.3);
Kop=sprandn(m,n,.1);
rhs1=rhs+sparse( Kop*ones(n,1)   );   % add and subtract at end
rhs=full(rhs);
tt=clock;
test=Kop'*sparse(rhs);
multtime=etime(clock,tt)
damp=0.;    % for +damp*norm(x)^2 in lss
atol=1000*eps;
btol=1000*eps;
conlim=1/atol;
itnlim=4*tn1;
msglvl=0;
aprod='my_aprod';
tt=clock;
[dir,se,itn,anorm,acond,rnorm,arnorm,xnorm,inform, ...
   ifail]=f04qaf(n,rhs1,aprod,damp, ...
             atol,btol,conlim,itnlim,msglvl);
dir=dir-ones(n,1);
n1s=[n1s n1];
ett=etime(clock,tt);
enorm=normest((Kop'*Kop*dir-Kop'*rhs));
%svdKop=svd(Kop);
disp(['dim and error and elapsed time is: ' ...
     num2str(n1),' ',num2str(enorm),'  ',num2str(ett)]);
etts=[etts ett];
multtimes=[multtimes multtime];
enorms=[enorms enorm];
end
plot(n1s,etts);
