LFSR in Sage

We may generate a keystream using the lfsr_sequence function, which takes as input the key (characteristic polynomial coefficients) and the initial state as vectors over a finite field, and the length of the keystream to be generated.

In [1]:
F = GF(2)
one = F.one()
key = [x*one for x in [1, 1, 1, 1]]
fill = [x*one for x in [0, 0, 0, 1]]
In [2]:
lfsr_sequence(key, fill, 20)
Out[2]:
[0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1]

Let us determine the periods of the LFSR.

In [3]:
def period(s, n):
    try:
        return next(i for i in range(1, len(s)) if s[i:i+n] == s[:n])
    except StopIteration:
        return None

def periods(key, lfsr=lfsr_sequence):
    F = key[0].parent()
    q = len(F)
    n = len(key)
    seen = []
    periods = {}
    l = q^n
    for fill in F^n:
        s = ''.join(str(x) for x in fill)
        if any(s in ss for ss in seen):
            continue
        seq = list(lfsr(list(key), list(fill), l+n))
        p = period(seq, n)
        periods[tuple(fill)] = p
        l -= p
        seen.append(''.join(str(x) for x in seq))
    return periods
In [4]:
periods(key)
Out[4]:
{(0, 0, 0, 0): 1, (1, 0, 0, 0): 5, (0, 1, 0, 0): 5, (1, 1, 1, 0): 5}

We observe that all nonzero states give a sequence of period 5. Let us verify that the connection polynomial $C(z)$ is irreducible and divides $z^5 - 1$.

In [5]:
P.<z> = PolynomialRing(F)
C = -1 + sum(a*z^i for i, a in enumerate(reversed(key), 1))
In [6]:
C.factor()
Out[6]:
z^4 + z^3 + z^2 + z + 1
In [7]:
(z^5 - 1) / C
Out[7]:
z + 1

Let us now consider a LFSR with the same recursion, only over a field of size 5.

In [8]:
F = GF(5)
one = F.one()
key = [x*one for x in [1, 1, 1, 1]]
In [9]:
periods(key)
Out[9]:
{(0, 0, 0, 0): 1, (1, 0, 0, 0): 312, (3, 1, 0, 0): 312}

We observe that all nonzero states give a sequence of period 312. Let us verify that the connection polynomial $C(z)$ is irreducible and divides $z^{312} - 1$.

In [10]:
P.<z> = PolynomialRing(F)
C = -1 + sum(a*z^i for i, a in enumerate(reversed(key), 1))
In [11]:
C.factor()
Out[11]:
z^4 + z^3 + z^2 + z + 4
In [12]:
(z^312 - 1) / C
Out[12]:
z^308 + 4*z^307 + 2*z^304 + 2*z^303 + z^302 + 4*z^300 + 2*z^299 + 4*z^297 + 3*z^296 + 3*z^294 + 3*z^293 + 2*z^292 + 2*z^291 + z^290 + 3*z^289 + z^288 + 2*z^287 + 4*z^284 + 3*z^283 + 3*z^282 + 3*z^280 + 2*z^279 + 3*z^278 + 2*z^277 + z^276 + z^275 + 4*z^274 + z^273 + z^271 + 2*z^270 + 3*z^269 + 4*z^268 + 2*z^267 + 3*z^266 + 4*z^265 + 4*z^262 + z^260 + 3*z^258 + z^257 + 2*z^256 + 4*z^255 + z^254 + 4*z^253 + 3*z^252 + z^251 + 3*z^250 + 2*z^249 + 2*z^248 + 4*z^247 + z^245 + 2*z^244 + z^243 + z^242 + 2*z^241 + 3*z^240 + z^238 + 3*z^237 + 4*z^236 + 2*z^235 + 2*z^234 + 2*z^230 + 3*z^229 + 4*z^226 + 4*z^225 + 2*z^224 + 3*z^222 + 4*z^221 + 3*z^219 + z^218 + z^216 + z^215 + 4*z^214 + 4*z^213 + 2*z^212 + z^211 + 2*z^210 + 4*z^209 + 3*z^206 + z^205 + z^204 + z^202 + 4*z^201 + z^200 + 4*z^199 + 2*z^198 + 2*z^197 + 3*z^196 + 2*z^195 + 2*z^193 + 4*z^192 + z^191 + 3*z^190 + 4*z^189 + z^188 + 3*z^187 + 3*z^184 + 2*z^182 + z^180 + 2*z^179 + 4*z^178 + 3*z^177 + 2*z^176 + 3*z^175 + z^174 + 2*z^173 + z^172 + 4*z^171 + 4*z^170 + 3*z^169 + 2*z^167 + 4*z^166 + 2*z^165 + 2*z^164 + 4*z^163 + z^162 + 2*z^160 + z^159 + 3*z^158 + 4*z^157 + 4*z^156 + 4*z^152 + z^151 + 3*z^148 + 3*z^147 + 4*z^146 + z^144 + 3*z^143 + z^141 + 2*z^140 + 2*z^138 + 2*z^137 + 3*z^136 + 3*z^135 + 4*z^134 + 2*z^133 + 4*z^132 + 3*z^131 + z^128 + 2*z^127 + 2*z^126 + 2*z^124 + 3*z^123 + 2*z^122 + 3*z^121 + 4*z^120 + 4*z^119 + z^118 + 4*z^117 + 4*z^115 + 3*z^114 + 2*z^113 + z^112 + 3*z^111 + 2*z^110 + z^109 + z^106 + 4*z^104 + 2*z^102 + 4*z^101 + 3*z^100 + z^99 + 4*z^98 + z^97 + 2*z^96 + 4*z^95 + 2*z^94 + 3*z^93 + 3*z^92 + z^91 + 4*z^89 + 3*z^88 + 4*z^87 + 4*z^86 + 3*z^85 + 2*z^84 + 4*z^82 + 2*z^81 + z^80 + 3*z^79 + 3*z^78 + 3*z^74 + 2*z^73 + z^70 + z^69 + 3*z^68 + 2*z^66 + z^65 + 2*z^63 + 4*z^62 + 4*z^60 + 4*z^59 + z^58 + z^57 + 3*z^56 + 4*z^55 + 3*z^54 + z^53 + 2*z^50 + 4*z^49 + 4*z^48 + 4*z^46 + z^45 + 4*z^44 + z^43 + 3*z^42 + 3*z^41 + 2*z^40 + 3*z^39 + 3*z^37 + z^36 + 4*z^35 + 2*z^34 + z^33 + 4*z^32 + 2*z^31 + 2*z^28 + 3*z^26 + 4*z^24 + 3*z^23 + z^22 + 2*z^21 + 3*z^20 + 2*z^19 + 4*z^18 + 3*z^17 + 4*z^16 + z^15 + z^14 + 2*z^13 + 3*z^11 + z^10 + 3*z^9 + 3*z^8 + z^7 + 4*z^6 + 3*z^4 + 4*z^3 + 2*z^2 + z + 1

lfsr_sequence only works with finite fields. We provide the function lfsr_ring which also lets us use ring elements.

In [13]:
def lfsr_ring(key, fill, n=-1):
    m = len(key)
    assert m == len(fill)
    i = 0
    for x in fill:
        yield x
        i += 1
        if i == n:
            return
    while i != n:
        x = sum([key[j] * fill[j] for j in range(m)])
        fill = fill[1:] + [x]
        yield x
        i += 1
In [14]:
R = Integers(10)
one = R.one()
key = [x*one for x in [1, 1, 1, 1]]
fill = [x*one for x in [0, 0, 0, 1]]
In [15]:
gen = lfsr_ring(key, fill)
seq = [next(gen) for _ in range(10000)]
print(seq[:100])
[0, 0, 0, 1, 1, 2, 4, 8, 5, 9, 6, 8, 8, 1, 3, 0, 2, 6, 1, 9, 8, 4, 2, 3, 7, 6, 8, 4, 5, 3, 0, 2, 0, 5, 7, 4, 6, 2, 9, 1, 8, 0, 8, 7, 3, 8, 6, 4, 1, 9, 0, 4, 4, 7, 5, 0, 6, 8, 9, 3, 6, 6, 4, 9, 5, 4, 2, 0, 1, 7, 0, 8, 6, 1, 5, 0, 2, 8, 5, 5, 0, 8, 8, 1, 7, 4, 0, 2, 3, 9, 4, 8, 4, 5, 1, 8, 8, 2, 9, 7]

Let us determine the period of a sequence.

In [16]:
period(seq, 4)
Out[16]:
1560
In [17]:
1560 / 5
Out[17]:
312

We see that the period of the sequence is the least common multiple of the periods appearing with the fields of sizes 2 and 5, which are prime factors of 10. In fact, every period will be obtained in this manner.

In [18]:
periods(key, lfsr_ring)
Out[18]:
{(0, 0, 0, 0): 1,
 (1, 0, 0, 0): 1560,
 (2, 0, 0, 0): 312,
 (5, 0, 0, 0): 5,
 (0, 1, 0, 0): 1560,
 (3, 1, 0, 0): 1560,
 (8, 1, 0, 0): 1560,
 (6, 2, 0, 0): 312,
 (0, 5, 0, 0): 5,
 (1, 1, 1, 0): 1560,
 (3, 1, 1, 0): 1560,
 (5, 5, 5, 0): 5}