% solve the H-equation (Kelley book, p87)
%      via Broyden's method
clear; close all

c = 1 - 1e-4;
N = input('(return for default) N = ');
if isempty(N), N = 2500; end
fprintf('H-equation: c = %g, N = %i\n\n',c,N)

x0 = ones(N,1);
tol = 1.e-12; 
maxiter = 100;
Lim = 4;
prt = 1;

t0 = tic;
[x,nrmF,iter] = yzBroyden(@yzHeqn,x0,tol,maxiter,prt,c,N);
fprintf('yzBroyden using yzHeqn: mean(x) = %.8f\n',mean(x))
fprintf('nrmF %g, iter %i, time %g\n\n',nrmF(end),iter,toc(t0))
h = semilogy(1:iter,nrmF,'-o'); set(h,'linewidth',2);
title('Convergence History'); 
xlabel('Iter'); ylabel('F-norm'); shg

if exist('myBroyden','file')
    t0 = tic;
    [x,nrmF,iter] = myBroyden(@yzHeqn,x0,tol,maxiter,prt,c,N);
    fprintf('myBroyden using myHeqn: mean(x) = %.8f\n',mean(x))
    fprintf('nrmF %g, iter %i, time %g\n\n',nrmF(end),iter,toc(t0))
    hold on, h = semilogy(1:iter,nrmF,'-o'); hold off
    set(h,'linewidth',2); legend('yzBroyden','myBroyden'); shg
end