% gradient descent quadratic bowl test
clc %clear console and return cursor to home
format long

%actual function
t = linspace(0,10,100);
N0 = 10;
R0 = 0.75;
Nd_smth = N0*exp(R0*t); 
M = length(Nd_smth);

%generate vals from random distr with mean = mn and stdv = stdv
rng('default'); %fix seed for random # generator to matlab startup seed
mn = 0; %mean
stdv = 1000;
r_pb = mn + stdv.*randn(1,M); 
Nd = r_pb+Nd_smth; %data
plot(t,Nd,'bx')

N = @(N0,R0,ti) N0.*exp(R0*ti); 
F = @(N0,R0,tm) 0.5*sum( (N(N0,R0,tm)-Nd).^2 );


F_N0 = @(N0,R0,tm) sum( (N(N0,R0,tm)-Nd).*exp(R0*tm) );
F_R0 = @(N0,R0,tm) sum( (N(N0,R0,tm)-Nd).*(N0*tm.*exp(R0*tm)) );

 % initial guess
 N0 = 100; R0 = 0.01; lambda = 1;
 %N0 = 100; R0 = 0.01; lambda = 50000;
 %N0 = 60; R0 = 0.51; lambda = 50000;
 %N0 = 27; R0 = 0.61;  lambda = 5000;
 %N0 = 15.8; R0 = 0.688;  lambda = 500;
 %N0 = 15.15; R0 = 0.7037; lambda = 500;
%  N0 = 13; R0 = 0.72; lambda = 1;
 
 
 %
 delta = 1e-9; %timestep or learning rate
 lambda = 10000;

 N_max = 60000; Tol = 1e-3;
 plot(t,Nd,'bx',t,N(N0,R0,t),'r-')
 pause;
 Out_max = 10;
 
for j = 1:Out_max
    for i = 1:N_max
        N0_old = N0;
        N0 = N0 - delta*(F_N0(N0,R0,t));
        R0 = R0 - (1/lambda)*delta*(F_R0(N0_old,R0,t));

        %plot curve every 5000 iterations
        if rem(i,5000) == 0
            N0
            R0
            i
            pause(0.1)
            plot(t,Nd,'bx',t,N(N0,R0,t),'r-')
        end        
    end
    %abs(F(N0,R0,t)) %output energy at current interation
    %pause   
        
        j
%         lambda = lambda/2;
%         if lambda < 1
%             lambda = 1;
%         end
        sprintf('lambda = %f', lambda)
        re = abs(F_obj(14.85,0.706,t,Nd)-F_obj(13,0.72,t,Nd))/abs(F_obj(13,0.72,t,Nd));
        sprintf('relative error = %f', re)
        %pause(2)
    
end
disp(N0)
disp(R0)

plot(t,Nd,'bx',t,N(N0,R0,t),'r-')

re = abs(F_obj(14.85,0.706,t,Nd)-F_obj(13,0.72,t,Nd))/abs(F_obj(13,0.72,t,Nd));
sprintf('relative error = %f', re)