Small Fibs Eventually Become Large Fibs – Another exercise in scaling

著者 Andy Campbell, May 31, 2016

2 ビュー (過去 30 日間) | 0 いいね | 6 コメント

OK, I had way too much fun writing that last blog post exploring algorithmic scaling. I want to explore this a bit more with a some "real" code, and what is more real than a function to calculate the lovely fibonacci number, the darling of interview questions, algorithms classes, and our own previous blogs. This will be fun to explore not only to cement the method we are using to explore scalability, but also because it is fun to see how fast and elegantly the fibonacci function can be written in MATLAB.

First, let's look at the most obvious solution based on the definition directly:

function f = fibonacci1(n)

if n < 1
    f = 0;
elseif n < 3
    f = 1;
else
    f = fibonacci1(n-1) + fibonacci1(n-2);
end

end

This is very nice and readable and reads exactly like the definition of fibonacci numbers, but we will see that unfortunately it is not very efficient at all. This is not efficient because the recursive approach needs to continually recalculate lower values of $n$ .

Let's create a parameterized test to invoke this algorithm. As I mentioned, this solution explodes pretty quickly so I am just going to take 10 linearly spaced points from 0 to 25. Also note that the calculation is fast for small n so we need to do this in a loop to get valid measurements.

classdef FirstFibonacciTest < matlab.perftest.TestCase
    
    
    properties(TestParameter)
        size = num2cell(round(linspace(0,25,10)));
    end
    
    
    methods(Test)
        function testFibonacci(~, size)
            for idx = 1:1000
                fibonacci1(size);
            end
        end
    end
end

How does it do?

recursiveResult = runperf('FirstFibonacciTest');

findMins = @(resultArray) arrayfun(@(result) min(result.Samples.MeasuredTime), resultArray);

ax = axes('XScale','log','YScale','log');
ylabel('Execution time (s)');
xlabel('Problem Size');
hold on;
grid on;

problemSize = linspace(0,25,10);
plot(ax, problemSize, findMins(recursiveResult), 'LineWidth', 2)
leg = legend({'Recursive'});

Running FirstFibonacciTest
..........
..........
..........
..........
..........
..........
..........
..........
.
Done FirstFibonacciTest
__________

$\mathcal{O}(2^n)$. Exponential. Problem sizes larger than 25 quickly become unproductive.

Well, what's the next approach? A common solution to calculate fibonacci numbers is to recognize that we don't need to recalculate $n-1$ and $n-2$ for every single $n$. Rather, we can just store these values into an array and calculate the number iteratively. Note this is actually a simple form of dynamic programming as it holds on to previous solutions of the problem in memory and it uses these prior solutions to calculate subsequent values. In essence we get some pretty substantial efficency gains because we leverage memory to our advantage. How does this look? We can just store the $nth$ fibonacci number along with those for $n+1$ and $n+2$ in a three element array and work through each of the fibonacci numbers iteratively continually updating our 3 element vector as we go.

function f = fibonacci2(n)

values = [0 1 1];
for q = 1:n
    values(1:2) = values(2:3);
    values(3) = sum(values(1:2));
end
f = values(1);

end

To observe how it scales I am going to modify our test to pass the algorithm in as another test parameter so we can reuse the test we have already written and ensure consistency of our measurements. We don't want one algorithm to be measured using a loop of 1000 and another accidently measured using a loop of 2000.

The test looks the same with the addition of the algorithm test parameter which contains function handles to our different fibonacci implementations.

classdef SecondFibonacciTest < matlab.perftest.TestCase
    
    
    properties(TestParameter)
        algorithm = struct(...
            'Recursive', @fibonacci1, ...
            'Iterative', @fibonacci2);
            
        size = num2cell(round(linspace(0,25,10)));
    end
    
    
    methods(Test)
        function testFibonacci(~, algorithm, size)
            for idx = 1:1000
                algorithm(size);
            end
        end
    end
end

Note that I have used the struct syntax to define the parameter which allows us to give nice names like "Recursive" and "Iterative" to each of the parameter values. This also allows us to select just these parameters and run the exact subset of interest, in this case, the iterative algorithm.

iterativeResult = runperf('SecondFibonacciTest', 'ParameterName', 'Iterative');
plot(ax, problemSize, findMins(iterativeResult), 'LineWidth', 2)
leg = legend([leg.String, 'Iterative']);

Running SecondFibonacciTest
..........
..........
..........
..........
..........
..........
..........
..........
...
Done SecondFibonacciTest
__________

Alright we now see the iterative approaches scales linearly ($\mathcal{O}(n)$)) which us much nicer, but as it turns out we can do better with a closed form solution using Binet's formula:

$$f_n = \frac{(1 + \sqrt{5})^n - (1 - \sqrt{5})^n}{2^n \sqrt{5}}$$

Here is how we code that up in MATLAB:

function f = fibonacci3(n)
sqrt5 = sqrt(5);
num = (1 + sqrt5)^n - (1 - sqrt5)^n; 
denom = (2^n)*sqrt5;
f = num/denom;
end

...and now we can simply add this new algorithm as new parameter to our existing test:

properties(TestParameter)
    size = num2cell(round(linspace(0,25,10)));
    algorithm = struct(...
        'Recursive', @fibonacci1, ...
        'Iterative', @fibonacci2, ...
        'ClosedForm', @fibonacci3);
end

and runperf it!

closedFormResult = runperf('FibonacciTest', 'ParameterName', 'ClosedForm');
plot(ax, problemSize, findMins(closedFormResult), 'LineWidth', 2)
leg = legend([leg.String, 'Closed Form']);

Running FibonacciTest
..........
..........
..........
..........
..........
..........
..........
..........
.
Done FibonacciTest
__________

That is looking very nice! In fact it demonstrates $\mathcal{O}(1)$ constant behavior. Can't get much better than that from a complexity perspective. However, there is one more solution I'd like to show using matrix exponentials. Hmmmm, I seem to remember hearing somewhere that MATLAB is pretty good with matrices. Does that sound familiar? Well lucky us, because the solution is quite trivial with the power of MATLAB. This is because the $nth$ fibonacci number $f_n$ is governed by the following equation:

$$ \left[ \begin{array}{cc} 1 & 1 \\ 1 & 0 \end{array} \right]^n = \left[ \begin{array}{cc} f_{n+1} & f_n \\ f_n & f_{n-1} \end{array} \right] $$

Yes, that's right, the solution in MATLAB is just that simple.

function f = fibonacci4(n)

F = [1 1;1 0]^n;
f = F(1,2);

end

How does this scale? Add the test and run it.

properties(TestParameter)
    size = num2cell(round(linspace(0,25,10)));
    algorithm = struct(...
        'Recursive', @fibonacci1, ...
        'Iterative', @fibonacci2, ...
        'ClosedForm', @fibonacci3, ...
        'MatrixExponential', @fibonacci4);
end

matrixResult = runperf('FibonacciTest', 'ParameterName', 'MatrixExponential');
plot(ax, problemSize, findMins(matrixResult), 'LineWidth', 2)
leg = legend([leg.String, 'Matrix Exponential']);

Running FibonacciTest
..........
..........
..........
..........
..........
..........
..........
..........
..........
..........
Done FibonacciTest
__________

Look at that! The closed form and the matrix solutions both show very scalable behavior, but the matrix form is not quite as good, showing logarithmic ($\mathcal{O}(\log n)$) scaling instead of constant and the constant factor is greater for the matrix case. However, there is another benefit enjoyed by the matrix exponential solution but not by the closed form solution. To show this, I am going to add a quick spot check to this test to confirm that we are actually getting the right answer. Note for this post I am just going to add another test method which validates a single size, but this could be done in the existing test method across all sizes if desired.

classdef VerifiedFibonacciTest < matlab.perftest.TestCase
    
    
    properties(TestParameter)
        size = num2cell(round(linspace(0,25,10)));
        algorithm = struct(...
            'Recursive', @fibonacci1, ...
            'Iterative', @fibonacci2, ...
            'ClosedForm', @fibonacci3, ...
            'MatrixExponential', @fibonacci4);           
    end
    
    
    methods(Test)
        function testFibonacci(~, algorithm, size)
            for idx = 1:1000
                algorithm(size);
            end
        end
        function testCorrectAnswer(testCase, algorithm)
           testCase.verifyEqual(algorithm(16), 987, ...
               'The 16th fibonacci number should be 987'); 
        end
    end
end

Let's run it with runtests to confirm correctness.

result = runtests('VerifiedFibonacciTest')

Running VerifiedFibonacciTest
..........
..........
..........
..........
..
================================================================================
Verification failed in VerifiedFibonacciTest/testCorrectAnswer(algorithm=ClosedForm).

    ----------------
    Test Diagnostic:
    ----------------
    The 16th fibonacci number should be 987

    ---------------------
    Framework Diagnostic:
    ---------------------
    verifyEqual failed.
    --> The values are not equal using "isequaln".
    --> Failure table:
                Actual    Expected           Error               RelativeError    
                ______    ________    ____________________    ____________________
            
                987       987         4.54747350886464e-13    4.60736930989325e-16
    
    Actual double:
             9.870000000000005e+02
    Expected double:
           987

    ------------------
    Stack Information:
    ------------------
    In /mathworks/inside/files/dev/tools/tia/mdd/performanceScaling/VerifiedFibonacciTest.m (VerifiedFibonacciTest.testCorrectAnswer) at 21
================================================================================
..
Done VerifiedFibonacciTest
__________

Failure Summary:

     Name                                                           Failed  Incomplete  Reason(s)
    ============================================================================================================
     VerifiedFibonacciTest/testCorrectAnswer(algorithm=ClosedForm)    X                 Failed by verification.
    

result = 

  1x44 TestResult array with properties:

    Name
    Passed
    Failed
    Incomplete
    Duration
    Details

Totals:
   43 Passed, 1 Failed, 0 Incomplete.
   7.3549 seconds testing time.

...and there is the problem with the closed form solution. It is subject to floating point precision errors. Of course we can simply make this test pass by leveraging a tolerance in verifyEqual

testCase.verifyEqual(algorithm(16), 987, 'RelTol', 1e-14...
    'The 16th fibonacci number should be 987');

but the fact remains that using the closed form solution we compute something that is close to the $nth$ fibonacci number, but the matrix solution computes it exactly. Which one to use depends on your application, but given the constant factor is so small and logarithmic behavior approaches constant for large problem sizes I think I would tend to choose the exact solution. How about you?

Published with MATLAB® R2016a