Categories::

• Root (97)
  • cryptography (25)
  • debian (13)
  • dev (18)
  • openbsd (2)
  • personal (10)
  • sage (25)

People::

Projects::

Bubble::

Gröbner Bases over $ZZ$

Vanilla Sage does not compute Gröbner bases over $ZZ$ by any definition. However, this feature has been requested several times. The earliest account I could find quickly is this post by Joe Wetherell. Below is a list of options for Gröbner bases over $ZZ$ in Sage.

1. Singular’s upcoming release will feature Gröbner bases over rings. In fact, the feature is present in the current Singular release but not enabled by default. An SPKG with that functionality enabled can be found here. ring r = integers,(x,y),lp; declares a ring over the integers where some things work and some things don’t. Note that this ring declaration is not final, i.e. the name integers may change. Also, this SPKG has issues and crashes on me for some operations. We’re working on tracking that issue down.
2. Macaulay2 has support for Gröbner bases over rings and a decent Sage interface supporting that functionality. Macaulay2 1.1 is available as an experimental SPKG. First, one needs to install boehm_gc-7.1.p0 and gdbm-1.8.3. Then, since the version in experimental didn’t compile for me, try my new SPKG.
3. If you are lucky enough to have Magma installed, I.groebner_basis(“magma:GroebnerBasis”) does the job. If you don’t have Magma installed try magma_free.
4. Ginv (also available as a optional package) also supports Gröbner bases over $ZZ$. However, for the example Joe gave in his e-mail it crashes on me and I’ve contacted upstream about it.
5. Last and least: I have a toy implementation of the $d$-Gröbner basis algorithm from the Becker-Weispfenning. Don’t hold your breath, it is dead slow.

Hopefully, due to the upcoming Singular release the situation will improve soon and we’ll finally have Gröbner bases over $ZZ$ in Sage.

posted at: 14:39 :: permanent link

Parallel Matrix Elimination

I released a new version of M4RI today which contains a parallel implementation for matrix elimination. Below I reproduce some timings for this code to give a rough idea of the performance of this code.

The examples hfe25_5, hfe30_5 and hfe35_5 from the M4RI website take 1.44, 9.29 and 51.56 seconds respectively.

Note that this is work in progress and that the algorithm still has worse complexity than the one implemented in Magma. Also note that the speed-up is far from linear and that the speed-up decreases with the size. This is probably because each thread falls out of L2 more often and the threads clog each other.

posted at: 22:56 :: permanent link

GGC 4.3 and -O3

I recently upgraded an Opteron server to Debian/Lenny to get GCC 4.3 for OpenMP reasons. It turns out that my code, namely matrix multiplication as implemented in the M4RI library, ran much slower than when compiled with GCC 4.1. For instance, to multiply two $20,000 \times 20,000$ random matrices took 18.38 seconds with GCC 4.1 but 21.00 seconds with GCC 4.3.1 and to multiply two $32,000 \times 32,000$ random matrices took 70.24 seconds with GCC 4.1 but 80.00 second with GCC 4.3.1. Eventually, I checked the highlevel changelog and found: “The -ftree-vectorize option is now on by default under -O3. In order to generate code for a SIMD extension, it has to be enabled as well: use -maltivec for PowerPC platforms and -msse/-msse2 for i?86 and x86_64.” However, we don’t use SSE2 on the Opteron since it is slower than the standard instruction set for this application. Passing -no-tree-vectorize to the compiler fixed the problem. However, to my surprise -O2 didn’t come with a speed penalty either, so I settled for this. The final timings on my Opteron server are:

I suppose the moral of the story is: -O3 isn’t necessarily better than -O2 just because 3>2.

posted at: 16:00 :: permanent link

ISSAC

I’m going to this year’s International Symposium on Symbolic and Algebraic Computation (ISSAC) starting on Sunday. William is giving a plenary lecture on Sage titled “Can There be a Viable Free Open Source Alternative to Magma, Maple, Mathematica, and Matlab?” a.k.a. “Computer algebra community meet Sage; Sage meet the computer algebra community”.

posted at: 14:23 :: permanent link

Scapy and Sage

“Scapy is a powerful interactive packet manipulation program. It is able to forge or decode packets of a wide number of protocols, send them on the wire, capture them, match requests and replies, and much more. It can easily handle most classical tasks like scanning, tracerouting, probing, unit tests, attacks or network discovery (it can replace hping, 85% of nmap, arpspoof, arp-sk, arping, tcpdump, tethereal, p0f, etc.). It also performs very well at a lot of other specific tasks that most other tools can’t handle, like sending invalid frames, injecting your own 802.11 frames, combining technics (VLAN hopping+ARP cache poisoning, VOIP decoding on WEP encrypted channel, …)”

At the end of the day Scapy is one (one!) Python file so it couldn’t be easier to use it from within Sage. As an example let’s assume we have sniffed an SSH connection establishment including a Diffie-Hellmann Group Exchange as described in RFC 4419. Scapy can do live packet capture and injection but that would require root privileges, so I’m working with a pcap file in this example:

from scapy import rdpcap, TCP, IP

SSH2_MSG_KEX_DH_GEX_GROUP = 31

# read packets
packets = [p[IP] for p in rdpcap("/home/malb/example.pcap") \
               if p[TCP] and len(p[TCP]) > 32]

# find correct package & payload
for packet in packets:
    try:
        pl = [ord(e) for e in packet[TCP].payload.load]
        if pl[5] == SSH2_MSG_KEX_DH_GEX_GROUP:
            break
    except AttributeError:
        pass

def get_uint(pl, length):
    # this is not as generic as it should be since it doesn't work
    # with negative numbers
    value = ZZ(0)
    for i in range(length):
        value += pl[i] * 2**(8*(length - i - 1))
    return value, pl[length:]

packet_length, pl = get_uint(pl, 4)
padlen, pl = get_uint(pl, 1)
packet_type, pl = get_uint(pl, 1)
assert(packet_type == SSH2_MSG_KEX_DH_GEX_GROUP)

# p
p_length, pl = get_uint(pl, 4)
p, pl = get_uint(pl, p_length)

# g
g_length, pl = get_uint(pl, 4)
g, pl = get_uint(pl, g_length)

assert(len(pl) == padlen)
assert(p.is_prime())

Zp = GF(p)
g = Zp(g)
e = g**ZZ.random_element(0,p)

e.log(g) # yeah, right ;-)

Happy hacking.

posted at: 17:34 :: permanent link

XOR for Fun and Profit

I just gave a talk on linear algebra over GF(2), optimisation techniques and applications to algebraic cryptanalysis. Slides are available online.

posted at: 22:02 :: permanent link

libM4RI in Debian Unstable

malb@XXX:~$ apt-cache search m4ri  
libm4ri-dev - Method of the Four Russians library, development files
libm4ri0 - Method of the Four Russians library, shared library

Big thanks to Tim for making that happen!

posted at: 06:37 :: permanent link

QOTD

“Mathematics is the art of reducing problems to linear algebra.” Bill Hart was introduced to this definition when he started his undergrad.

posted at: 02:12 :: permanent link

Fraction Free Gauss-Jordan Errata

I’m at Sage’s dev1 right now and so I have the pleasure of meeting Arne Storjohann. In his thesis he presented a fraction free asymptotically fast matrix elimination algorithm which unfortunately has some typos in it. Below I replicate the correct algorithm he explained/provided me yesterday:

def GaussJordan(A, k=-1, d0=None):
  if d0 == None:
    d0 = A.base_ring()(1)
  n = A.nrows()
  m = A.ncols()
  I = MatrixSpace(A.base_ring(),A.nrows(),A.nrows())(1)
  for i in xrange(k+1,n):
    if any(A[i,j] for j in xrange(m)):
      break
  else:
      U,P,r,h,d = d0*I, I, 0, n-k, d0
      return (U,P,r,h,d)

  if m == 1:
    i = min([i for i in xrange(k+1,n) if A[i,0] != 0])
    P = copy(I)
    P.swap_rows(i,k+1)
    r,h,d = 1, n-i, (P*A)[k+1,0]
    U = d*I
    for j in range(n):
      U[j,k+1] = -(P*A)[j,0]
    U[k+1,k+1] = d0
  else:
     m1,m2 = m//2, m-m//2
     A1 = A.matrix_from_columns(range(m1))
     B = A.matrix_from_columns(range(m1,m))
     U1, P1, r1, h1, d1 = GaussJordan(A1, k, d0)
     A2 = d0**(-1)*U1*P1*B
     U2, P2, r2, h2, d = GaussJordan(A2, k+r1, d1)
     U = d1**(-1) * U2*(P2*(U1 - d1*I) + d1*I)
     P,r,h,d = P2*P1, r1+r2, min(h1,h2),d
  return U,P,r,h,d

Note that this is not how one would actually implement this algorithm in practice: it is pseudo-code that happens to run in Sage. For a practical implementation check the IML library.

posted at: 17:56 :: permanent link

AES Equation Systems

From time to time either I or Carlos receive requests for AES and/or BES equation systems. This post is an attempt to score high enough on Google so that others can find out about the AES equation system generator that is in Sage by default. For a while now Sage shippes a generator for AES and BES equation systems and their small scale variants (called “SR”). The generator supports both $GF(2)$ and $GF(2^e)$ for $e \in \{4,8\}$.

# Help for generator.
sage: mq.SR?
# We construct SR(1,1,1,4) over GF(2^4).
sage: sr = mq.SR(1,1,1,4)
SR(1,1,1,4)
# The constructor may fail due to zero inversions.
sage: F,s = sr.polynomial_system()
---------------------------------------------------------------------------
<type 'exceptions.ZeroDivisionError'>     Traceback (most recent call last)
...
<type 'exceptions.ZeroDivisionError'>: A zero inversion occurred during an 
encryption or key schedule.
# So we try again.
sage: F,s = sr.polynomial_system(); F
Polynomial System with 40 Polynomials in 20 Variables
# The object F is a polynomial system and the object s a solution
# dictionary. Help about F can be found via tab completion.
sage: F.<tab>

# F can be exported to Magma.
sage: sage: F._magma_()
Ideal of Polynomial ring of rank 20 over GF(2^4)
Graded Reverse Lexicographical Order
Variables: k100, k101, k102, k103, x100, x101, x102, x103, w100, w101, w102, 
w103, s000, s001, s002, s003, k000, k001, k002, k003
Basis:
[
w100 + k000 + $.1^4,
w101 + k001 + $.1^8,
w102 + k002 + $.1,
w103 + k003 + $.1^2,
k000^2 + k001,
....
....
....

# F can be exported to Singular.
sage: F._singular_()
w100+k000+(a+1),
w101+k001+(a^2+1),
w102+k002+(a),
w103+k003+(a^2),
k000^2+k001,
k001^2+k002,
k002^2+k003,
...

# Or we can use those systems transparently in the background.
sage: F.groebner_basis() # Singular in the background
[k002 + (a^3 + 1)*k003 + (a^2), k001 + (a^3 + a^2)*k003 + (a^3), k000 + 
(a^2)*k003 + (a^3 + a^2), s003 + (a^3 + a)*k003 + (a^3 + a^2 + a), s002 + 
(a^3 + a^2 + a)*k003 + (a^2), s001 + (a^3 + a^2 + a + 1)*k003 + (a + 1), s000 
+ (a^2 + a)*k003 + 1, w103 + k003 + (a^2), w102 + (a^3 + 1)*k003 + (a^2 + a), 
w101 + (a^3 + a^2)*k003 + (a^3 + a^2 + 1), w100 + (a^2)*k003 + (a^3 + a^2 + a 
+ 1), x103 + (a^3 + a)*k003, x102 + (a^3 + a^2 + a)*k003 + (a^3 + 1), x101 + 
(a^3 + a^2 + a + 1)*k003 + (a^3 + a), x100 + (a^2 + a)*k003 + (a^3 + a), k103 
+ (a^3 + a + 1)*k003 + 1, k102 + (a^2 + a + 1)*k003 + (a^3 + a^2), k101 + (a 
+ 1)*k003 + (a + 1), k100 + (a)*k003 + (a^2 + a + 1), k003^2 + (a^2)*k003 + 
(a^3 + a^2)]
sage: F.groebner_basis(algorithm='magma:GroebnerBasis') # Magma
[k003^2 + (a^2)*k003 + (a^3 + a^2), k100 + (a)*k003 + (a^2 + a + 1), k101 + (a 
+ 1)*k003 + (a + 1), k102 + (a^2 + a + 1)*k003 + (a^3 + a^2), k103 + (a^3 + a 
+ 1)*k003 + 1, x100 + (a^2 + a)*k003 + (a^3 + a), x101 + (a^3 + a^2 + a + 
1)*k003 + (a^3 + a), x102 + (a^3 + a^2 + a)*k003 + (a^3 + 1), x103 + (a^3 + 
a)*k003, w100 + (a^2)*k003 + (a^3 + a^2 + a + 1), w101 + (a^3 + a^2)*k003 + 
(a^3 + a^2 + 1), w102 + (a^3 + 1)*k003 + (a^2 + a), w103 + k003 + (a^2), s000 
+ (a^2 + a)*k003 + 1, s001 + (a^3 + a^2 + a + 1)*k003 + (a + 1), s002 + (a^3 
+ a^2 + a)*k003 + (a^2), s003 + (a^3 + a)*k003 + (a^3 + a^2 + a), k000 + 
(a^2)*k003 + (a^3 + a^2), k001 + (a^3 + a^2)*k003 + (a^3), k002 + (a^3 + 
1)*k003 + (a^2)]

# We can also construct equation systems over GF(2).
sage: sr = mq.SR(2,1,1,4,gf2=True)
sage: F,s = sr.polynomial_system()

# For those we can use PolyBoRi to compute the Groebner basis.
sage: R= F.ring()
sage: B = BooleanPolynomialRing(R.ngens(), R.variable_names(), order="lex")
sage: F2 = B.ideal([B(f) for f in F]) 
sage: F2.groebner_basis()
[k200 + k001 + k003 + 1, k201 + k001, k202 + 1, k203 + k000, x200 + k003, x201 
+ k000 + k001, x202 + k000 + k001 + k003, x203 + k000 + k003, w200 + k000 + 
k003 + 1, w201 + k001 + k003 + 1, w202 + k001 + 1, w203, s100 + k003 + 1, 
s101 + k000 + k001, s102 + k000 + k001 + k003, s103 + k000 + k003, k100 + 1, 
k101 + k001 + k003, k102 + k000 + 1, k103 + k003, x100 + k001 + k003 + 1, 
x101 + k000 + k001 + k003 + 1, x102 + k001 + 1, x103 + k000, w100 + k000 + 1, 
w101 + k001 + 1, w102 + k003, w103 + k003 + 1, s000 + k000 + k001 + k003, 
s001 + k000 + k001, s002 + k000 + k003, s003 + k001 + 1, k000*k001 + k000 + 
k001 + 1, k000*k003 + k000 + k003 + 1, k001*k003 + k001, k002 + k003]

Happy attacking.

posted at: 12:20 :: permanent link

M4RI Website

I finally put together the website for the M4RI library. For those who don’t know M4RI:

“M4RI is a library for fast arithmetic with dense matrices over $\mathbb{F}_2$. It was started by Gregory Bard and is now maintained by Martin Albrecht and Gregory Bard. The name M4RI comes from the first implemented algorithm: The “Method of the Four Russians” inversion algorithm published by Gregory Bard. This algorithm in turn is named after the “Method of the Four Russians” multiplication algorithm which is probably better referred to as Kronrod’s method. M4RI is used by the Sage mathematics software and the PolyBoRi library. M4RI is available under the General Public License Version 2 or later (GPLv2+).

Features of the M4RI library include:

• basic arithmetic with dense matrices over $\mathbb{F}_2$ (addition, equality testing, stacking, augmenting, sub-matrices, randomisation)
• asymptotically fast $O(n^{log_27})$ matrix multiplication via the “Method of the Four Russians” (M4RM) & Strassen-Winograd algorithm,
• asymptotically fast $O(n^{3}/log_2(n))$ row echelon form computation and matrix inversion via the “Method of the Four Russians” (M4RI), and
• support for the x86/x86_64 SSE2 instruction set where available.
• support for Linux and OS X (GCC), support for Solaris (Sun Studio Express) and support for Windows (Visual Studio 2008 Express).”

Performance-wise it is doing okay but not great. On Intel’s Core2Duo it seems to compare favourably to Magma 2.13. Though, I don’t have access to Magma 2.14 yet which improves dense linear algebra over $\mathbb{F}_2$. However, on AMD’s Opteron it is way behind Magma 2.13. This is possibly due to the 1MB L2 cache of the Opteron vs. 4MB L2 cache of the Core2Duo.

posted at: 12:26 :: permanent link

Algebraic Techniques in Differential Cryptanalysis

Finally, our paper on algebraic techniques in differential cryptanalysis is available as pre-print. From the abstract:

“In this paper we propose a new cryptanalytic method against block ciphers, which combines both algebraic and statistical techniques. More specifically, we show how to use algebraic relations arising from differential characteristics to speed up and improve key-recovery differential attacks against block ciphers in some situations. To illustrate the new technique, we apply it to reduced round versions of the cipher PRESENT, an ultra lightweight block cipher proposed at CHES 2007, particularly suitable for deployment in RFID tags.”

Since the results of that paper rely on experimental data, we also publish the source code used to execute our experiments. I’m going to present this paper at SCC 2008, whichs schedule is now available, btw.

posted at: 13:22 :: permanent link

A Cryptographic Tour and Todo List of Sage

Yesterday someone showed up on [sage-devel] and wrote: “I have been developing software and doing research in the areas of: mathematics, cryptography algorithms, encryption, and would like to contribute my time and effort to the Sage project. I would like any of you to get me started in the right direction, any info would be appreciated.”

This is the edited/polished version of my reply. I am posting it here in case anyone else wonders how to contribute to Sage for cryptographic research.

• David Kohel wrote an introductionary book to cryptography. He uses Sage in the book and wrote a fair amount of code to make that happen. The relevant module is sage.crypto. For example, it implements a linear feedback shift register:

  sage: FF = FiniteField(2)
  sage: P.<x> = PolynomialRing(FF)
  sage: E = LFSRCryptosystem(FF); E
  LFSR cryptosystem over Finite Field of size 2
  
  sage: IS = [ FF(a) for a in [0,1,1,1,0,1,1] ]
  sage: g = x^7 + x + 1
  sage: e = E((g,IS))
  sage: B = BinaryStrings()
  sage: m = B.encoding("THECATINTHEHAT")
  sage: e(m)
  0010001101111010111010101010001100000000110100010101011100001011110010010000011111100100100011001101101000001111
  

  It seems the code needs more documentation and also some areas are not implemented yet, e.g. block ciphers.
• Sage ships PyCrypto which implements many standard cryptographic algorithms. The docstring level documentation is horrible:

  sage: import Crypto.Cipher.IDEA
  sage: Crypto.Cipher.IDEA?
       x.__init__(...) initializes x; see x.__class__.__doc__ for signature
  
  

  It is not really meant for research/education/playing around but for production code but maybe something could be done to have easier access to it from within Sage. Here is an example how to use it:

  sage: from Crypto.Hash import MD5
  sage: m = MD5.new()
  sage: m.update('abc')
  sage: m.digest()
  '\x90\x01P\x98<\xd2O\xb0\xd6\x96?}(\xe1\x7fr'
  sage: m.hexdigest()
  '900150983cd24fb0d6963f7d28e17f72'
  

• Finite fields are basic building blocks in cryptography. Sage uses several finite field implementations for prime fields and extension fields of various sizes. The FiniteField_ext_pari implementation for finite extension fields of order $\ge 2^{16}$ should be replaced by two implementations using NTL’s ZZ_pE and lzz_pE depending on the size of the characteristic. This should be relatively straight-forward because there is an implementation for characteristic 2 using NTL’s GF2E already. To get a feeling about the possible speed improvements:

  sage: k.<a> = GF(next_prime(2^65)^27)
  sage: e = a^30
  sage: f = a^40
  sage: %timeit e*f
  1000 loops, best of 3: 557μs per loop
  
  sage: c = ntl.ZZ_pEContext(ntl.ZZ_pX(list(k.polynomial()),k.characteristic()))
  sage: e = c.ZZ_pE(list(e.polynomial()))
  sage: f = c.ZZ_pE(list(f.polynomial()))
  sage: %timeit e*f
  10000 loops, best of 3: 154μs per loop
  

• Sage isn’t exactly kicking ass when it comes to elliptic and hyperelliptic curves over finite fields. As these are quite important in asymmetric cryptography it might be worth looking into this. John Cremona added to this: “That is fair. Apart from an SEA point-counting implementation — only over prime fields — the rest (for elliptic curves) is definitely only designed to work at sub-crypto field sizes.”
• algebraic techniques received some attention for the cryptanalysis of symmetric cryptographic primitives recently. In these attacks the cryptanalyst expresses the cipher as a large set of multivariate polynomial equations and attempts to solve the system. The most common case over $\mathbb{F}_2$ is handled by PolyBoRi. This library is the backbone of BooleanPolynomialRing and friends. This class needs testing, documentation, extension and bugfixes. Basically someone should sit down and add all the methods of MPolynomial[Ring]_libsingular to BooleanPolynomial[Ring] which make sense, add a ton of doctests and test the hell out of the library to make sure no SIGSEGVs surprise the user.

  sage: F,s = sr.polynomial_system()
  sage: R = F.ring()
  sage: B = BooleanPolynomialRing(R.ngens(),R.variable_names(),R.term_order())
  sage: F = [B(f) for f in F if B(f) != 0]
  sage: F = mq.MPolynomialSystem(B,F); F
  Polynomial System with 68 Polynomials in 36 Variables
  sage: gb = F.groebner_basis()
  
  sage: gb[-1]
  k003
  sage: s[R("k003")]
  0
  

• The module sage.crypto.mq is also relevant for algebraic cryptanalysis in symmetric cryptography. It implements
  • small scale AES equation system generators over $\mathbb{F}_2$ and $\mathbb{F}_{2^n}$
  • a class to represent multivariate polynomial systems
  • an S-box class to analyse … well … S-boxes.

  sage: S = mq.SBox(7, 6, 0, 4, 2, 5, 1, 3); S
  (7, 6, 0, 4, 2, 5, 1, 3)
  
  sage: S.polynomials()
  
  [x0*x2 + x1 + y1 + 1,
   x0*x1 + x1 + x2 + y0 + y1 + y2 + 1,
   x0*y1 + x0 + x2 + y0 + y2,
   x0*y0 + x0*y2 + x1 + x2 + y0 + y1 + y2 + 1,
   x1*x2 + x0 + x1 + x2 + y2 + 1,
   x0*y0 + x1*y0 + x0 + x2 + y1 + y2,
   x0*y0 + x1*y1 + x1 + y1 + 1,
   x1*y2 + x1 + x2 + y0 + y1 + y2 + 1,
   x0*y0 + x2*y0 + x1 + x2 + y1 + 1,
   x2*y1 + x0 + y1 + y2,
   x2*y2 + x1 + y1 + 1,
   y0*y1 + x0 + x2 + y0 + y1 + y2,
   y0*y2 + x1 + x2 + y0 + y1 + 1,
   y1*y2 + x2 + y0]
  
  sage: S.difference_distribution_matrix()
  
  [8 0 0 0 0 0 0 0]
  [0 2 2 0 2 0 0 2]
  [0 0 2 2 0 0 2 2]
  [0 2 0 2 2 0 2 0]
  [0 2 0 2 0 2 0 2]
  [0 0 2 2 2 2 0 0]
  [0 2 2 0 0 2 2 0]
  [0 0 0 0 2 2 2 2]
  

• Univariate polynomials over $\mathbb{F}_2$ are still implemented via NTL’s ZZ_pX rather than GF2X.

  sage: k = GF(2)
  sage: P.<x> = PolynomialRing(k)
  sage: f = P([k.random_element() for _ in range(10000)])
  sage: %timeit f**2
  100 loops, best of 3: 11.1 ms per loop
  sage: f = ntl.GF2X(f.list())
  sage: %timeit f**2
  100000 loops, best of 3: 2.02μs per loop
  

  Furthermore there is gf2x a drop-in replacement library for NTL’s GF2X which is expected to be five times faster than NTL. Though, a formal vote is needed to get it into Sage.
• At the end of the day everything boils down to linear algebra. So if one improves that, everybody wins. Sparse linear algebra over $\mathbb{F}_p$ is still too slow (Ralf-Phillip Weinmann did some work here wrapping code from eclib), there is no special implementation for sparse linear algebra over $\mathbb{F}_2$ (both blackbox and reduced echelon forms), dense linear algebra over $\mathbb{F}_2$ lacks Strassen multiplication/reduction and dense linear algebra over $\mathbb{F}_{2^n}$ should probably get a specialised implementation. Note that up until a few thousand rows/columns that Sage’s dense linear algebra over $\mathbb{F}_2$ is actually pretty fast.

  sage: A = random_matrix(GF(2),7000,7000)
  sage: time B = A.echelon_form()
  CPU times: user 2.91 s, sys: 0.04 s, total: 2.96 s
  Wall time: 3.01
  
  sage: AM = A._magma_()
  sage: t = magma.cputime()
  sage: BM = AM.EchelonForm()
  sage: magma.cputime(t)
  2.6800000000000002
  

I hope this list isn’t totally useless.

posted at: 14:06 :: permanent link

Yet Another Talk on Sage

I gave a brief talk yesterday at the Open Knowledge Conference (OKCon) here in London. The slides were also discussed on [sage-devel] last week. I have to admit that I underappreciated David Joyner’s comments about the expected audience. My impression is that the majority of the audience couldn’t care less about the actual mathematics implemented in Sage. I suppose we still made a good impression but I had to skip most of the examples I care about due to time constraints and preceived lack of interest. After the talks I had some neat discussions with other participants, e.g. Gaël Varoquaux from the MayaVi2 project.

posted at: 18:44 :: permanent link

Seriously Cool

This gets even me excited about physics:

via: http://www.badscience.net/?p=622

posted at: 23:19 :: permanent link

Plotting Timing Experiments

Like any other person I regulary need to run experiments to check how fast or slow a particular algorithm/implementation is for a given problem. The natural choice is to plot the data. This way you at least get some more or less pretty picture out of the tendious experience of having to wait for the experiment to finish. I used to write crappy code to generate these pictures myself and I could not convince myself to remember the appropriate commands for matplotlib and R. Today I sat down and learned the five lines of code necessary to have decent plots for my experiments. I’m putting examples here for no good reason except maybe to show off Sage’s new HNF code which I use as a showcase.

First lets compare how long it takes to compute the Hermite Normal Form for a given random $n \times n$ matrix with (possibly negative) integer entries of size bounded absolute by $2^{16}$.

n = 10
b = 16
st =[]
mt = []
x = [20*i for i in range(n)]
for i in range(n):
  A = random_matrix(ZZ,20*i,20*i, x=-2**b, y=2**b)
  t = cputime()
  E = A.echelon_form()
  st.append(cputime(t))

  AM = A._magma_()
  t = magma.cputime()
  EM = AM.EchelonForm()
  mt.append(magma.cputime(t))

import pylab

pylab.clf() # clear the figure first

pylab.figure(1)

# plot some data and add a legend
pylab.plot(x,st,label="Sage") 
pylab.plot(x,mt,label="Magma") 
pylab.legend() # print the legend

pylab.title("HNF for Random Matrices with $%d$-bit Integer Entries: Sage vs. Magma"%b)
pylab.ylabel("execution time $t$") # label the axes
pylab.xlabel("n for n x n matrix")

pylab.savefig('foo.png',dpi=72) # fire!

Now lets use R to see how the runtimes vary for random $160 \times 160$ matrices with (possible negative) integer entries bounded absolute by $2^{10}$.

b = 10
st = []
for i in range(500):
  A = random_matrix(ZZ,160,160, x=-2**b, y=2**b)
  t = cputime()
  E = A.echelon_form()
  st.append(cputime(t))

from rpy import r

r.png('histogram.png',width=640,height=480)
r.hist(st,r.seq(1.2,3.7,0.02),main="SAGE HNF Histogram",col="lightblue", prob=True, xlab="seconds")
r.lines(r.density(st,bw=0.05),col="black")
r.rug(st)
r.dev_off()

Neat, isn’t it? Btw. Pygments is also neat, thanks rpw.

posted at: 23:08 :: permanent link

“Algebraic Techniques in Differential Cryptanalysis”

That is the title of the talk I am going to give tomorrow at the ISG Research Seminar:

Abstract: We propose a new cryptanalytic method against block ciphers, which combines both algebraic and statistical techniques. More specifically, we show how to use algebraic relations arising from differential characteristics to speed up and improve key-recovery differential attacks against block ciphers in some situations. To illustrate the new technique, we apply it to reduced round versions of the cipher PRESENT, an ultra lightweight block cipher proposed at CHES 2007, particularly suitable for deployment in RFID tags.

posted at: 14:28 :: permanent link

Impressions from FSE 2008

If you don’t get it, don’t worry, it is not really funny.

posted at: 14:17 :: permanent link

CFP: First International Conference on Symbolic Computation and Cryptography

I was asked to forward this announcement a while ago and finally get around to do so:

“SCC 2008 is the first of a new series of conferences where research and development in symbolic computation and cryptography may be presented and discussed. It is organized in response to the growing interest in applying and developing methods, techniques, and software tools of symbolic computation for cryptography. The use of lattice reduction algorithms in cryptology and the application of Groebner bases in the context of algebraic attacks are typical examples of explored applications.

SCC 2008 aims at providing an interactive forum for interested researchers to exchange ideas and views, to present research results and progress, and to learn and discuss recent developments and emerging problems on

• the design, modeling, and analysis of cryptographic systems and protocols for which symbolic computation may be used or needed, and
• the design, implementation, and analysis of algorithms and software
• tools of symbolic computation that may have potential applications in cryptography.

TOPICS

Specific topics for SCC 2008 include, but are not limited to:

• Multivariate cryptography, braid group cryptography, noncommutative cryptography, and quantum cryptography
• Code-based, factorization-based, and lattice-based cryptography
• Algebraic attacks for block ciphers, stream ciphers, and hash functions
• Design and analysis of algebraic, elliptic, and embedded cryptographic systems and protocols
• Groebner basis techniques in cryptology, algebraic number theory, and coding theory
• Triangular sets and new techniques for solving algebraic systems over finite fields
• Algorithms and software for symbolic computation in cryptography

INVITED SPEAKERS

Bruno Buchberger (Johannes Kepler Universitaet Linz, Austria)
Arjen K. Lenstra (Ecole Polytechnique Federale de Lausanne, Switzerland) (to be confirmed)
Adi Shamir (Weizmann Institute of Science, Israel)
Xiaoyun Wang (Tsinghua University and Shandong University, China)
…”

posted at: 13:58 :: permanent link

Les Trophées du Libre 2007

Sage is among the finalists of this year’s “free software awards” competition in the science category. The other two finalists in that category are Giac/XCas (slides, session) and Getfem++. I am representing Sage in 25 minutes and I uploaded my slides and the demo worksheet (PDF).

posted at: 11:53 :: permanent link

$GF(2^n)$ arithmetic speed

Since version 2.8.10 Sage’s finite extension fields of characteristic 2 and degree $\ge 16$ are implemented via NTL’s GF2E rather than Pari. For some more or less random reason I timed how fast multiplying two random elements is now.

The red line show the time it takes Magma 2.13-5 to multiply two random elements a million times for a given degree $n$. The green line shows the same calculation using Sage 2.8.12 with the default modulus and a Python loop. The blue line uses a Cython loop (== C loop) and the function good_modulus (see below) to generate a “good” modulus. The default modulus used by Sage is either the conway polynomial or - if we don’t know the conway polynomial - a random irreducible polynomial. I took the idea of using a “good” modulus from Michael Scott’s slides for his talk at the SPEED workshop. My attempt is not as sophisticated as his but naively searches for trinomials and pentanomials with low degree terms.

def good_modulus(n):
  P = GF(2)['x']
  x = P.gen()
  for a in xrange(1,n):
    f = x**n + x**a + 1
    if f.is_irreducible():
      return f
  for N in range(0,n,10):
    for a in xrange(1,N+1):
      for b in xrange(a+1,N+1):
        for c in xrange(b+1,N+1):
          f = x**n + x**c + x**b + x**a + 1
          if f.is_irreducible():
            return f
  # fall back to default if nothing was found
  return GF(2**n,'a').polynomial()

Some comments:

1. Up until $2^{15}$ we use Zech logarithms as they are implemented in Givaro. Magma uses Zech logarithms up to $2^{20}$ and we should do the same. If we use a Cython loop (i.e. remove the overhead of the loop) Sage’s arithmetic is as fast as Magma’s.
2. I don’t know why there is that peak around $n=2$ for Magma. Bug? My bad?
3. Magma scales quite nicely wordwise, as you would expect.
4. Surprisingly enough we beat Magma starting at $2^{100}$ up until at least $2^{128}$ using the “good” moduli.
5. What is going on with NTL between $2^{16}$ and $2^{64}$?

It seems we should internally - at least for large degrees - represent elements w.r.t. to a “good” modulus even if we know the conway polynomial.

posted at: 18:36 :: permanent link

Yet Another Talk on Sage

May I point the reader’s attention to the slides of my most recent talk about Sage for the ISG Student Seminar.

posted at: 19:40 :: permanent link

More Pictures/SAGE Days 5

I’ve uploaded my pictures from the “Tools for Cryptanalysis 2007” and “SAGE Days 5” workshops to flickr. At SD5 I

• worked on LLL (see wiki, trac, and update),
• found out that SAGE’s Strassen Echelonizer doesn’t require submatrices to be non-singular,
• generated kinda handy graphs of SAGE’s inheritance tree,
• and gave a talk about the status of commutative algebra in SAGE.

The SAGE inheritance tree in 3D:

posted at: 15:53 :: permanent link

Small Scale AES MQ Generator

I’ve just filled a SAGE trac ticket which has my implementation of a small scale AES polynomial system generator attached. It supports all SR sizes, systems over $GF(2)$ and $GF(2^e)$, and $SR*$ both as specified and in AES mode, i.e. such that $SR*(10,4,4,8)$ is AES. The patch also contains a base class called MPolynomialSystemGenerator which is supposed to be a common plattform for similar constructors for other ciphers. I would appreciate feedback on that. The generator spits out an instance of MPolynomialSystem which I’ve written of before: “I am planing to provide a class for SAGE for polynomial systems (over $GF(2)$ and $GF(2^e)$) and would appreciate input what kind of constructions should be possible with this class. As a working title I call it MPolynomialSystem (done) and it will provide methods to construct stuff like coefficient matrices, ideals etc.(done) This by itself is not that exciting, but I am planing to subclass this to allow more fine grade access to e.g. block cipher rounds. (done) Also, MPolynomialSystem classes over $GF(2^e)$ should be convertible to ideal bases over GF(2) (done) and so on.”

A little example session:

sage: sr = mq.SR(1,1,1,4,gf2=True)
sage: sr
SR(1,1,1,4)
sage: P = sr.random_state_array(); P
[a^2 + a + 1]
sage: K = sr.random_state_array(); K
[a + 1]
sage: sr(P,K)
[a^2 + a + 1]
sage: F,s = sr.polynomial_system(P,K)
sage: F
Polynomial System with 56 Polynomials in 20 Variables
sage: type(F)
<class 'sage.crypto.mq.mpolynomialsystem.MPolynomialSystem_gf2'>
sage: s
{k003: 1, k002: 1, k001: 0, k000: 0}
sage: gb = F.groebner_basis()

sage: sr = mq.SR(1,1,1,4,gf2=False)
sage: F,s = sr.polynomial_system(P,K)
sage: type(F)
<class 'sage.crypto.mq.mpolynomialsystem.MPolynomialSystem_gf2e'>
sage: gb = F.groebner_basis()
sage: F2 = F.change_ring(GF(2))
sage: F2
Polynomial System with 240 Polynomials in 80 Variables
sage: type(F2)
<class 'sage.crypto.mq.mpolynomialsystem.MPolynomialSystem_gf2'>

Comments? Feedback? Rants?

posted at: 19:16 :: permanent link

Upcoming Workshops/ Travel Plans

I am going to visit:

• ECrypt’s Tools for Cryptanalysis Workshop, Institute of Mathematics (IMPAN), Krakow (Poland), September 24 - 25, 2007.
• SAGE Days 5, Clay Math Institute, Cambridge (USA), September 30 - October 3, 2007.
• SAGE Days 6, Heilbr