% USEMIN_SMART  Script to test convergence of MINIMAL.  Starts each
% new iteration with an interpolated version of last solution.
% ELB 4/24/05

% grid sizes to test
JJ=[10 20 40 80 120];
TOL=1e-6;

%  g and f for manufactured example
g=inline('4*x.*y','x','y');
f=inline('128*x.*y.*(1+16*(x.^2+y.^2)).^(-1.5)','x','y');

N=length(JJ);  err=zeros(1,N);  time=zeros(1,N);
for k=1:N
   J=JJ(k);  
   if k==1
      xx=0:1/J:1;  [XX,YY]=ndgrid(xx,xx); 
      uex=g(XX,YY);
   else
      xx=0:1/J:1;  [XXnew,YYnew]=ndgrid(xx,xx);
      uex=g(XXnew,YYnew);
      U0=interpn(XX,YY,U,XXnew,YYnew,'cubic');  % interpolate previous U
      XX=XXnew;  YY=YYnew;
   end
   disp(' ')
   disp(['START minimal with J = ' num2str(J) ',  TOL = ' num2str(TOL) ':'])
   if k==1
      tic, U=minimal(J,g,TOL,f); time(k)=toc;
   else
      tic, U=minimal(J,g,TOL,f,U0); time(k)=toc;
   end
   err(k)=max(max(abs(U-uex)));
   disp(['DONE.  Results:  err = ' num2str(err(k)) ',  time = ' num2str(time(k))])
end

figure
subplot(1,2,1), loglog(1./JJ,err,'o')
xlabel('\Delta x'),  ylabel('error'),  set(gca,'XTick',sort(1./JJ))
axis([0.8*min(1./JJ) 1.2*max(1./JJ) 0.8*min(err) 1.2*max(err)])
subplot(1,2,2), loglog(JJ,time,'*')
xlabel J, ylabel('execution time'),  set(gca,'XTick',JJ)
axis([0.8*min(JJ) 1.2*max(JJ) 0.8*min(time) 1.2*max(time)])