Wednesday 30 April 2014

Fast Matrix Exponentiation


The code typed-in during the tutorial-session is as follows:-

    //Buttons: Number Exp, Fib Num Gen

    using System;
    using System.Collections.Generic;
    using System.ComponentModel;
    using System.Data;
    using System.Diagnostics;
    using System.Drawing;
    using System.Linq;
    using System.Text;
    using System.Threading.Tasks;
    using System.Windows.Forms;

    namespace Tutorial
    {
        public partial class Form1 : Form
        {
            public Form1()
            {
                InitializeComponent();
            }

            //Exponentiation by Squaring
            //4^16                =    (4^2)^8
            //16^8                =    (16^2)^4
            //256^4              =    (256^2)^2
            //65536^2          =    (65536^2)^1
            //4294967296    =    4294967296

            //4^15                        =    4 * (4^2)^7                   =   4 * 16^7
            //4 * 16^7                  =    4 * 16 * (16^2)^3         =   64 * 256^3
            //64 * 256^3              =    64 * 256 * (256^2)^1   =   16384 * 65536^1
            //16384 * 65536^1    =    16384 * 65536
            //1073741824

            private static Int64 NumberExponentiation(Int64 num, Int64 power)
            {
                if (power == 0)
                    return 1;

                if (power == 1)
                    return num;

                if (power % 2 == 0)
                    return NumberExponentiation(num * num, power / 2);
                else
                    return num * NumberExponentiation(num * num, (power - 1) / 2);
            }

            private void button1_Click(object sender, EventArgs e)
            {
                Int64 num = 11;
                Int64 power = 5;

                //11^5 = 161051.
                MessageBox.Show(NumberExponentiation(num, power).ToString());
            }

            private static Int64 TraditionalFibGen(Int64 termNum)
            {
                //Fibonacci Series:-
                //1, 1, 2, 3, 5, 8, 13, 21, 34...

                Int64 prev = 1;
                Int64 curr = 1;
                Int64 fibNum = 0;

                for (Int64 i = 3; i <= termNum; i++)
                {
                    fibNum = prev + curr;
                    prev = curr;
                    curr = fibNum;
                }

                if (fibNum == 0)
                    return 1;
                else
                    return fibNum;
            }

            private void button2_Click(object sender, EventArgs e)
            {
                Int64 termNum = 100000000;

                Stopwatch sw = Stopwatch.StartNew();
                Int64 tradFibNum = TraditionalFibGen(termNum);
                sw.Stop();

                string tradFibStatus = "Trad Fib Num: " + tradFibNum + "\nElapsed Ticks: " + sw.ElapsedTicks;
                //MessageBox.Show(tradFibStatus);

                //Fibonacci Series:-
                //1, 1, 2, 3, 5, 8, 13, 21, 34...

                //Matrix Method:-
                //The count of numbers that don't depend on any previous number(s) is 2.
                Int64[][] initialMatrix = new[] { new[] { 1L }, new[] { 1L } };

                //Matrix Multiplication:-
                //  _        _     _        _        _                            _       _          _
                // |   1  2   | * |   2  1   |  =  |   1*2+2*2    1*1+2*1   |  =  |   6   3   |
                // |_ 3  4 _|    |_ 2  1 _|      |_ 3*2+4*2    3*1+4*1 _|      |_ 14  7 _|

                //We need to find what transformation matrix when multiplied with the initial matrix yields the next set of Fib numbers.

                //axb * 2x1 = 2x1
                //b = 2
                //a = 2

                //The dimension of initialMatrix is 2x1.
                //The dimension of the matrix after matrix-multiplication must be 2x1.
                //Hence, the transformation matrix must have a 2x2 dimension.

                //  _        _     _    _        _    _
                // |   a  b   | * |   1   |  =  |   1   |
                // |_ c  d _|    |_ 1 _|      |_ 2 _|

                // Equations:-
                // a*1 + b*1 = 1
                // c*1 + d*1 = 2

                //Hence, a = 0; b = 1; c = 1; d = 1.
                Int64[][] transformation = new[] { new[] { 0L, 1L }, new[] { 1L, 1L } };

                //5th Fib Number Generation:-
                //  _        _                                  _    _      _    _
                // |   0  1   | raised to power 4 * |   1   | = |   5   |
                // |_ 1  1 _|                                 |_ 1 _|    |_ 8 _|

                //This is not element-wise powering. It's matrix-multiplication powering.

                //Detailed example for arriving at the 5th Fib number:-
                //  _        _     _        _        _                            _       _         _
                // |   0  1   | * |   0  1   |  =  |   0*0+1*1    0*1+1*1   |  =  |  1   1   |
                // |_ 1  1 _|    |_ 1  1 _|      |_ 1*0+1*1    1*1+1*1 _|      |_ 1   2 _|

                //  _      _       _        _       _                            _       _        _
                // |   1  1   | * |  0  1   |  =  |   1*0+1*1    1*1+1*1   |  =  |  1   2  |
                // |_ 1  2 _|    |  1  1 _|      |_ 1*0+2*1    1*1+2*1 _|      |_ 2   3_|

                //  _        _     _        _        _                            _       _          _
                // |   1  2   | * |   0  1   |  =  |   1*0+2*1    1*1+2*1   |  =  |   2   3   |
                // |_ 2  3 _|    |_ 1  1 _|      |_ 2*0+3*1    2*1+3*1 _|      |_ 3   5 _|

                //  _         _     _     _        _         _       _    _
                // |   2   3   | * |   1   |  =  |   2*1+3*1   |  =  |   5   |
                // |_ 3   5 _|    |_ 1 _|      |_ 3*1+5*1 _|      |_ 8 _|

                Int64 termNum2 = termNum - 1;

                sw = Stopwatch.StartNew();

                Int64[][] finalMatrix = MatrixMultiplication(MatrixExponentiation(transformation, termNum2), initialMatrix);
                sw.Stop();

                string fastFibStatus = "\nFast Fib Num: " + finalMatrix[0][0] + "\nElapsed Ticks: " + sw.ElapsedTicks;
                MessageBox.Show(tradFibStatus + fastFibStatus);
            }

            //  _        _      _        _       _                           _       _          _
            // |   1  2   | * |   2  1   |  =  |  1*2+2*2    1*1+2*1   |  =  |   6   3   |
            // |_ 3  4 _|    |_ 2  1 _|      |_3*2+4*2    3*1+4*1 _|      |_ 14  7 _|

            private static Int64[][] MatrixMultiplication(Int64[][] a, Int64[][] b)
            {
                int aRows = a.GetLength(0);
                int aCols = a[0].GetLength(0);
                int bRows = b.GetLength(0);
                int bCols = b[0].GetLength(0);

                Int64[][] c = new Int64[aCols][];

                for (int i = 0; i < aCols; i++)
                    c[i] = new Int64[bRows];

                for (int i = 0; i < aRows; i++)
                    for (int j = 0; j < bCols; j++)
                        for (int k = 0; k < aCols; k++)
                            c[i][j] = (c[i][j] + a[i][k] * b[k][j]);
               
                return c;

            }

            private static Int64[][] MatrixExponentiation(Int64[][] matrix, Int64 power)
            {
                if (power == 0)
                    return new[] { new[] { 0L, 1L }, new[] { 1L, 1L } };

                if (power == 1)
                    return matrix;

                if (power % 2 == 0)
                    return MatrixExponentiation(MatrixMultiplication(matrix, matrix), power / 2);
                else
                    return MatrixMultiplication(matrix, MatrixExponentiation(MatrixMultiplication(matrix, matrix), (power - 1) / 2));
            }

        }
    }

Thank you for reading this post and watching the videos. Please Subscribe, Comment, and Rate the channel of you liked the videos.

Goto C# Experiments to access more of such content! Thanks again!


No comments:

Post a Comment