Fibonacci LFSRs calculation optimisation

Question

Fibonacci LFSR is described on wiki, it's pretty simple.
I'd like to calucate the period of some Fibonacci's LFSR and use generated sequence for ciphering later.
Let's take and example from wiki:

x¹⁶ + x¹⁴ + x¹³ + x¹¹ + 1;

//code from wiki:
include <stdint.h>
uint16_t lfsr = 0xACE1u;
unsigned bit;
unsigned period = 0;
do {
  /* taps: 16 14 13 11; characteristic polynomial: x^16 + x^14 + x^13 + x^11 + 1 */
  bit  = ((lfsr >> 0) ^ (lfsr >> 2) ^ (lfsr >> 3) ^ (lfsr >> 5) ) & 1;
  lfsr =  (lfsr >> 1) | (bit << 15);
  ++period;
} while(lfsr != 0xACE1u);

My weakly try so far in php:

function getPeriod(){
    $polynoms = array(16, 14, 13, 11);
    $input = $polynoms[0] - 1;
    $n = sizeof($polynoms);
    for ($i = 1; $i < $n; $i++)
        $polynoms[$i] = $polynoms[0] - $polynoms[$i];
    $polynoms[0] = 0;
    //reversed polynoms == array(0, 2, 3, 5);
    $lfsr = 0x1; //begining state       
    $period = 0;
    //gmp -- php library for long numbers;
    $lfsr = gmp_init($lfsr, 16);
    do {            
        $bit = $lfsr; //bit = x^16 >> 0;

        for($i = 1; $i < $n; $i++) {
            //bit ^= lfsr >> 2 ^ lfst >> 3 ^ lfst >> 5;
            $bit = gmp_xor($bit, ( gmp_div_q($lfsr, gmp_pow(2, $polynoms[$i])) ));                  
        }
        //bit &= 1;
        $bit = gmp_and($bit, 1);            
        //lfsr = $lfsr >> 1 | $lsfr << (16 - 1);
        $lfsr = gmp_or( (gmp_div_q($lfsr, 2)), (gmp_mul($bit, gmp_pow(2, $input))) );

        $period++;
    } while (gmp_cmp($lfsr, 0x1) != 0);
    echo '<br />period = '.$period;
    //period == 65535 == 2^16 - 1; -- and that's correct;
            //                        I hope, at least;
    return $period;
}

Problem:
If i try to modulate work of i.e.

x³²¹ + x¹⁴ + x¹³ + x¹¹ + 1;

i got an error:"Fatal error: Maximum execution time of 30 seconds exceeded in /var/www/Dx02/test.php";

Can i somehow optimize (accelerate :) ) the calculation?

Any help is appreciated. Thank you and excuse me for my English.

Answer 1

To optimize this a bit we need to remember that PHP has great overhead on parsing code as it is not compiled, so we need to do as much work ourselves for it as we can. You should always profile your CPU/memory sensitive code with xdebug+KCachegrind(for example) to see where PHP wastes most of it's time. With your code only 12% is spent on gmp_* functions calculations, most of the time is spent for code parsing.

On my notebook(it is rather slow) my code runs 2.4 sec instead of 3.5 sec for your code, but for greater degrees the difference should be more noticeable (for example 19 power gives 19 vs 28 sec). It is not much, but it is some.

I left comments inside code, but if you have some questions - feel free to ask. I used function creation to replace that 'for($i = 1; $i < $n; $i++)' loop inside your main loop.

Also, I think you should change type of your $period variable to GMP(and $period++ to gmp_* function) as it can be greater then maximum integer on your system.

function getPeriod() {
    $polynoms = array(16, 14, 13, 11);


    $highest = $polynoms[0];
    $input = $highest - 1;

    //Delete first element of array - we don't need it anyway
    array_shift($polynoms);

    $polynoms_count = count($polynoms);

    //You always repeat gmp_pow(2, $input) and it's result is constant,
    //so better precalculate it once.
    $input_pow = gmp_pow(2, $input);

    //Start function creation.
    //If you don't use PHP accelerators, then shorter variable names
    //work slightly faster, so I replaced some of names
    //$perion->$r,$bit -> $b, $lfsr -> $l, $polynoms -> $p
    $function_str = '$r=0;';
    $function_str .= 'do{';

    //Now we need to get rid of your loop inside loop, we can generate
    //static functions chain to replace it.
    //Also, PHP parses all PHP tokens, even ';' and it takes some time,
    //So, we should write as much one-liners as we can.
    $function_str .= '$b=gmp_xor($b=$l';


    foreach ($polynoms AS $id => &$polynom) {
        //You always repeat gmp_pow(2, $polynoms[$i]) and it's result is constant,
        //so better precalculate it once.
        $polynom = gmp_pow(2, $highest - $polynom);

        //We create our functions chain here
        if ($id < $polynoms_count - 1) {
            $function_str.=',gmp_xor(gmp_div_q($l,  $p[' . $id . '])';
        } else {
            $function_str.=',gmp_div_q($l,  $p[' . $id . '])';
        }
    }

    //Close all brackets
    $function_str.=str_repeat(')', $polynoms_count);

    //I don't know how to optimize the following, so I left it without change
    $function_str.=';';
    $function_str.='$l = gmp_or((gmp_div_q($l, 2)), (gmp_mul(gmp_and($b, 1), $i_p)));';
    $function_str.='$r++;';
    $function_str.='} while (gmp_cmp($l, 0x1));';
    $function_str.='return $r;';

    //Now, create our funciton
    $function = create_function('$l,$p,$i_p', $function_str);


    //Set begining states
    $lfsr = 0x1;
    $lfsr = gmp_init($lfsr, 16);


    //Run function
    $period = $function($lfsr, $polynoms, $input_pow);

    //Use result
    echo '<br />period = ' . $period;
    return $period;
}

Answer 2

You simply can't do it this way with a polynomial like x^321 + ...;

If the polynomial is chosen well, you get a period length of 2^231 -1, and this is approximately 4.27 *10^96. If I'm not mistaken, this number is believed to exceed the number of atoms in the universe... (Strictly speaking, I'm referring to the posted C-code since I do not know php, but that certainly makes no difference.)

However, there is a mathematical method to calculate the length of the period without doing a brute-force attack. Unfortunately, this can't be explained in a few lines. If you have a solid background in math (especially calculations in finite fields), I'll be glad to look for a helpful reference for you.

EDIT: The first step in calculating the period of the LFSR obtained by using a polynomial p(x) is to obtain a factorization of p(x) mod 2, i.e. in GF(2). To do this, I recommend using software like Mathematica or Maple if available. You could also try the freely available Sage, see e.g. http://www.sagemath.org/doc/constructions/polynomials.html for usage details.

The period of p(x) is given by its order e, that means the smallest number such that p(x) divedes x^e+1. Unfortunately, I can't provide more information at the moment, it will take me several days to look for the lecture notes of a course I took several years ago...

A small example: p(x) = (x^5+x^4+1) = (x^3+x+1)*(x^2+x+1), the individual periods are 2^3-1=7 and 2^2-1=3, and since all polynomial factors are different, the period of p(x) is 3*7=21, which I also verified in C++.