function [x,e,xt,et] = BR_KB_Jacobian(KB,KBF,alpha,tol,max_iter,x0)
% BR_KB_Jacobian - gradient descent to solve PSAT
% On input:
%     KB (KB structure): knowledge base
%     KBF (string): KB function name
%     alpha (float): step size
%     tol (float): termination tolerance
%     max_iter (int): max number of iterations to try
%     x0 (vector): initial guess
% On output:
%     x (nx1 vector): solution
%     e (float): sentence error generated by solution
%     xt (kxn array): trace of descent
%     et (kx1 vector): trace of sentence error
% Call:
%     [x,e,xt,et] = BR_KB_Jacobian(KB1,'KB1i',0.001,0.001,1000,rand(20,1));
% Author:
%     T. Henderson
%     UU
%      Spring 2018
%

%figure(1);
%clf
%hold on

BETA = 0.8;

n = length(x0);
done = 0;
iter = 0;
x = x0;
xt = x0';
%et = feval(KBF,x0);
vt = feval(KBF,x0);
et = norm(vt);
e = et(end);
wb = waitbar(0,'Jacobian');
del_x = 0;
mu = 0.1;
while done==0
    alpha = BETA*alpha;
    waitbar(iter/max_iter);
    iter = iter + 1;
    J = zeros(n,1);
    for d = 1:n
        [iter d];
%        plot(iter,d,'k.');
%        drawnow
        xm = x;
        xm(d) = x(d) - alpha;
%        vm = feval(KBF,xm);
        vm = norm(feval(KBF,xm));
        xp = x;
        xp(d) = x(d) + alpha;
%        vp = feval(KBF,xp);
        vp = norm(feval(KBF,xp));
        J(d) = vm - vp;
    end
    if e>tol&norm(J)==0
        J = 0.5 - rand(n,1);
    end
    J = J/norm(J);
    del_x = mu*del_x + alpha*J;
    x = abs(x+del_x);
%    x = abs(x + alpha*J);
    x = max(0,min(x,1));
    xt = [xt; x'];
    et = [et;norm(feval(KBF,x))];
    e = et(end);
    if et(end)<tol||(iter==max_iter)
        done = 1;
    end
end
close(wb);
e = et(end);
