%% Introduction to Programming
% This part of the course provides a brief example-driven introduction to computer programming in Matlab.  We will focus on the computer programming concepts that are helpful for solving dynamic models.
% We will discuss two types of loops that are helpful for solving a dynamic model: "while" and "for" loops.
%
%% While Loops
% Suppose we needed an algorithm to solve the sqare root of two.  This is an irrational number, so it is not easy to calculate with a pencil and paper.  Instead, we rely on a so-called while loop to help us approximate it.
% The algorithm will proceed as follows: we know that 1^2 = 1 and 2^2 = 4.  This means that the square root of two must be bounded between 1 and 2.  One easy thing to do is to check the midpoint between 1 and 2 (i.e., 1.5) and square that to see if it's close to 2:

1.5^2

% We see that 1.5^2 = 2.25.  Not a bad guess.  Now we know that the square root of two must be bounded between 1 and 1.5.  Instead of calculating this manually over and over again, we'll just turn this job over to a while loop.
% A while loop will continue evaluating a set of lines until a particular condition is met.  In this case, it makes sense to end this cycle once the squared approximation is sufficiently close to two.  We'll call this sufficiency a tolerance, since an irrational number can only be approximated.  We'll call the difference between our approximation and our target value the error.   Let's give it a shot:

tolerance   = 0.000001 ;
upper_bound = 4        ;
lower_bound = 1        ;
error       = 10       ;    % We're just picking something large here to get the while loop going.

iteration   = 0        ;    % Let's keep track of the number of iteration it takes to get to our target.

while error > tolerance    % Keep iterating while this is true.  The loop will stop when this is not true any more.
    
    approximation = (upper_bound + lower_bound) / 2 ;    % Taking the midpoint, or "bisecting" the interval.
    
    error         = abs(approximation^2 - 2) ;    % Measuring how close we are to our target.
    
    % Now we test the approximation and update the interval bounds.
    
    if approximation^2 < 2    % In this case, our lower bound is too low.
        
        lower_bound = approximation ;    % Square root of two must be between the approximation and the upper bound.
        
    elseif approximation^2 >= 2
        
        upper_bound = approximation ;    % Square root of two must be between the approximation and the lower bound.
        
    end
    
    iteration = iteration + 1;    % Updating our iteration counter.
    
end

fprintf( 'The square root of two is approximately %0.5f.\n'  , approximation )
fprintf( 'We solved this in %d iterations.\n'                , iteration     )

%
% In the numerical analysis world, this method is known as bisection.  It can be used to find the root of a function, and sometimes it's used in numerical optimization to find where a function's derivative is zero.
% Notice that we introduced a new concept: the if/else statement.  We will use this frequently to navigate the computer through the code when two or more mutually exclusive and mutually exhaustive evaluations are available.
% Be careful with while loops - notice that if we had forgotten to evaluate the error, then this loop would have run infinitely.  Generally want to limit the number of evaluations of a while loop to ensure this doesn't happen.  This can be done by setting the maximum number of iterations and adding a condition to the while loop:

maximum_iterations = 10000;
while ( error > tolerance ) && ( iteration < maximum_iterations )
end

%%
%% For Loops
% The second type of loop we will discuss is called a for loop.  Generally, we'll use for loops when we have to do a specific number of calculations that can be mapped to the natural numbers.
% Suppose, for example, that we wanted to evaluate f(x) = x^2 at several points between -10 and 10.  To do that, we could choose, say 15 points, and space them linearly.  Then for each of these points\{x_1, x_2,...,x_{15}\}, we evaluate f(x_i).
% To evaluate f(x_i) at each of these 15 points, we could just use a for loop (i.e., "for" each of the 15 calculations):

% Assigning parameters and initalizing values
upper_bound =  10                                           ;
lower_bound = -10                                           ;
n_evals     =  15                                           ;
x_values    = linspace( lower_bound, upper_bound, n_evals ) ;
fx_values   = zeros( 1, n_evals )                           ;    % Preallocating the output vector to ensure Matlab won't have to do a costly real-time expansion.


for ie = 1:n_evals
    
    fx_values(ie) = x_values(ie)^2 ;
    
end

plot(x_values, fx_values)

% Generally, we would like to avoid for loops because they usually take up lots of space, and they can run more slowly than alternative operations.
% A solution to this is vectorization.  This invoves using linear algebra and element-wise operations to do the heavy lifting.
% For example, we could have simplified the previous code and gained efficiency by simply writing:

% Assigning parameters and initalizing values
upper_bound =  10                                           ;
lower_bound = -10                                           ;
n_evals     =  15                                           ;
x_values    = linspace( lower_bound, upper_bound, n_evals ) ;
% fx_values   = zeros( 1, n_evals )                           ;    % No need to preallocate here, because the values are assigned all at once.

fx_values = x_values.^2 ;    % Emphasis on the ".^" operation.  Generally, a "." before an operator will make it an element-wise operation.

plot(x_values, fx_values)

%
%