I computed M74207281

On 7th of January 2016, GIMPS found Mersenne prime number M74207281. GIMPS is a project combining efforts of a large crowd of volunteers contributing computing power and a university server combining the result – maintained by Curtis Cooper. If a number can be found satisfying the conditions of a Mersenne prime number, the number will be checked with primality tests; according to the article it took 31 days to do so.

M74207281 is specified 2274 207 281-1 and a friend of mine at GoGraz wondered how to compute this number in decimal representation. I suggested square and multiply, but wasn’t sure how long it might take.

So I implemented the exp-by-squaring-iterative algorithm on Wikipedia (didn’t thoroughly check whether it is the most efficient algorithm for this kind of problem):

package main

import (
	"log"
	"math/big"
	"os"
	"strconv"
)

var zero *big.Int
var one *big.Int
var two *big.Int
var hundred *big.Int

func init() {
	// initialize global values
	zero = big.NewInt(0)
	one = big.NewInt(1)
	two = big.NewInt(2)
	hundred = big.NewInt(100)
}

// ExpBySquaringIterative exponentiates by squaring meaning it computes x^n
// in logarithmic time in size of n
func ExpBySquaringIterative(x, n *big.Int) *big.Int {
	switch n.Cmp(zero) {
	case -1:
		x.Div(one, x)
		n.Neg(n)
	case 0:
		c := big.NewInt(1)
		return c
	}

	y := big.NewInt(1)
	for n.Cmp(one) > 0 {
		x2 := big.NewInt(0)
		x2.Set(x)

		nm := big.NewInt(0)
		if n.Bit(0) == 0 {
			x.Mul(x2, x2)
			n.Div(n, two)
		} else {
			y.Mul(x2, y)
			x.Mul(x2, x2)
			nm.Sub(n, one)
			n.Div(nm, two)
		}
	}
	return x.Mul(x, y)
}

func main() {
	if len(os.Args) == 1 {
		log.Println("usage: go run num_2n1.go <exponent>")
		log.Println("  M74207281 is given by exponent 74207281")
		os.Exit(1)
	}

	expstr := os.Args[1]
	exponent, err := strconv.Atoi(expstr)
	if err != nil {
		log.Fatal(err)
	}

	log.Println("TBD   ExpBySquaringIterative")
	result := ExpBySquaringIterative(big.NewInt(2), big.NewInt(int64(exponent)))
	result.Sub(result, one)
	log.Println("DONE  ExpBySquaringIterative")
	log.Println("TBD   computing decimal representation")
	log.Printf("2^%d - 1 = ", exponent)
	log.Println(result)
	log.Println("DONE  computing decimal representation")
}

And how long does it take to compute this number?

Computation of M74207281

Figure 1. Computation of M74207281

Was there anything interesting about it?

  • The decimal representation uses 22.4 MB memory.
  • It took 1 hour 27 minutes to compute the decimal representation. But it only takes 4 seconds to actually compute the number in binary. By far most of the hours are spent on computing the decimal string representation.
  • The builtin function Exp (you need to set m to nil) takes 1 hour 17 minutes; so it takes a little bit less time.
  • For reference, the first digits are 300376418084… and the last digits are …073391086436351 as can be seen on lcns2’s page.
I computed M74207281

When is Tutte’s theorem in K₅ violated?

Theorem

A graph, G = (V, E), has a perfect matching if and only if for every subset U of V, the subgraph induced by V − U has at most |U| connected components with an odd number of vertices.
Tutte’s theorem

Illustration
tutte_theorem_examples
Figure 1. Tutte theorem examples

Consider the two examples given in Figure 1. G is a bipartite graph and we consider two possible subsets U the induced subgraph of G-U is given on the right side.

  1. For U₁ the number of odd components in G-U is 1 and the size of U is given with 3. Because 1 ≤ 3, the Tutte theorem holds.
  2. For U₂ the number of odd components in G-U is also 1 and the size of U is also given with 3. Again, the Tutte theorem holds.

For any U, this condition will hold. Underneath the perfect matching is given (vertex 1 is associated to vertex 2, etc).

Problem setting

planar graph K5

Figure 2. Planar graph K5 and a reduced variant

Consider graph K₅. A perfect matching clearly cannot exist for K₅, because a perfect matching can only exist if the number of vertices is even (as every vertex is associated to another vertex). Because the Tutte theorem provides a necessary and sufficient condition for the existence of a perfect matching, it should yield

For some subset U of V, the subgraph induced by V − U
has less than |U| connected components with an odd number of vertices.

Question

What is U for K₅? The solution is equivalent for its incomplete variant, which was originally given to me.

Solution

I wasn’t able to come up with a solution. So I wrote a program to try all combinations:

#!/usr/bin/env python3

import itertools

def component_sizes(G):
    """component sizes in an undirected graph"""
    comps = list(set([v]) for v in G[0])
    for e in G[1]:
        s, t = e[0], e[1]
        idx_s, idx_t = -1, -1
        for i, comp in enumerate(comps):
            if s in comp:
                idx_s = i
            if t in comp:
                idx_t = i
        if idx_s == idx_t:
            continue
        else:
            for item in comps[idx_t]:
                comps[idx_s].add(item)
            comps.remove(comps[idx_t])
    return tuple(sorted(len(comp) for comp in comps))

assert component_sizes(([1], [])) == (1,)
assert component_sizes(([1, 2], [])) == (1, 1)
assert component_sizes(([1, 2], [(1, 2)])) == (2,)
assert component_sizes(([1, 2, 3], [(1, 2)])) == (1, 2)

def subgraph(G, excl_vs):
    """Subgraph induces by excluding vertices"""
    red_V = tuple(filter(lambda v: v not in excl_vs, G[0]))
    red_E = tuple(filter(lambda e: e[0] not in excl_vs and e[1] not in excl_vs, G[1]))
    return (red_V, red_E)

assert subgraph(([1], []), [1]) == (tuple(), tuple())
assert subgraph(([1, 2], [(1, 2)]), [1]) == ((2,), tuple())
assert subgraph(([1, 2, 3], [(1, 2), (2, 3), (3, 1)]), [1]) == ((2, 3), ((2, 3),))

if __name__ == '__main__':
    def excl_vs_selection(G):
        for size in range(len(G[0])):
            for sel in itertools.combinations(G[0], r=size):
                yield set(sel)

    K5 = ((1, 2, 3, 4, 5), list(itertools.combinations(range(1,6), r=2)))
    reduced_K5 = ((1, 2, 3, 4, 5), [(1, 2), (2, 3), (3, 4), (4, 5), (5, 1)])

    G = K5
    for excl_vs in excl_vs_selection(G):
        sub = subgraph(G, excl_vs)
        sizes = component_sizes(sub)
        odd_sized_comps = len(list(map(lambda v: v % 2 == 0, sizes)))
        if odd_sized_comps > len(excl_vs):
            print(excl_vs)
            #raise ValueError("Tutte theorem violated with exclusion vertices {}".format(excl_vs))

The program yields set(). Hence the only solution is the empty set.

Consider U = {}. Then the induced graph is K₅ itself. The number of connected components in K₅ is 1 (all pairs of vertices are reachable from each other) and this component has an odd number of vertices (namely 5). The size of U is 0. And 1 > 0, so the Tutte theorem is violated and no perfect matching exists indeed.

When is Tutte’s theorem in K₅ violated?

Make 24 out of 1, 3, 4 and 6

Problem

Give a mathematical expression which equals 24. This expression is only allowed to use the basic arithmetic operators +, -, * and /. The numbers 1, 3, 4 and 6 must occur exactly once. No magic or tricks!

Setting

I was given this problem last Friday by one colleague. It was ridiculously difficult to come up with a solution. As mathematician you try to understand the structure behind the problem. You look at the prime numbers involved. In the expression you have {3, 2, 2, 2, 3} and in 24 you have {2, 2, 2, 3}. So you need to get rid of one three. But how? None of us 6 guys came up with a solution in less than an hour. But we actually found one academic who found the solution after 2 hours (a bit of sleep was involved, because he was very tired).

Well, my mathematical brain failed, but my computer science brain screamed “brute force!”. So let’s do it. I wrote a python script, which solved a subproblem in 18 lines. This is a more complete version which solves the problem for any four numbers and any target value. This script prints possible correct solutions to stdout in less than 0.5 seconds:

import itertools

def parenthesized(a, b, c, d):
    """Generate all expressions where operators
    are represented with substitution characters
    """
    yield '((`{} $ <{}) @ (>{} ! _{}))'.format(a, b, c, d)
    yield '((`{} $ (<{} @ >{})) ! _{})'.format(a, b, c, d)
    yield '(`{} $ ((<{} @ >{}) ! _{}))'.format(a, b, c, d)
    yield '(`{} $ (<{} @ (>{} ! _{})))'.format(a, b, c, d)

def main():
    """Main routine generating math expressions equating 24
    using numbers {1, 3, 4, 6} exactly once and using operators
    -, +, / or *.
    """
    for arr in itertools.permutations([1, 3, 4, 6]):
        for expr in parenthesized(*arr):
            # we can omit '-' in these iterations, because all values
            # will also be considered with a sign
            for fo, so, to in itertools.product(['+', '/', '*'], repeat=3):
                for s1, s2, s3, s4 in itertools.product(['', '-'], repeat=4):
                    try:
                        # replace operators
                        e = expr.replace('$', fo).replace('@', so).replace('!', to)
                        # use either sign '' (hence +) or '-' (hence -) at every place
                        e = e.replace('`', s1).replace('<', s2).replace('>', s3).replace('_', s4)
                        result = eval(e)
                    except ZeroDivisionError:
                        continue
                    if result == 24:
                        yield e

if __name__ == '__main__':
    for solution_expression in main():
        print(solution_expression)
Spoiler

In the following I will give a list of hints. The further below you go, the closer you come to the solution (select text to display items). The solution itself is not given here:

  1. Most values we ended up with are too small for 24. So we conjectured that number 6 will occur as factor, which turned out to be true.
  2. Factor 6? So we look for factor 4 as well, right? Actually no. I mean yes, but with an indirect representation.
  3. Given {1, 4} on the left-hand side and {1, 3, 4} on the right-hand side. Find operators such that LHS and RHS equate.
  4. How about 1 / (1 / 4)? Two ones? I guess one of them has to be 6.
  5. Given 6 / (1 / 4)? How can we represent 1/4 using {1, 3, 4}?
Make 24 out of 1, 3, 4 and 6

Running Donald Knuth’s CWEB programs

I just realized that no-one documented how to quickly run Donald Knuth’s programs he publishes in his WEB system.

Suppose you want to compile a program like commafree; a program he wrote between September and December 2015 for his Annual Christmas tree lecture in 2015 about Universal Commafree Codes.

  1. Download the webfile which comprises documentation and source code.

    wget http://www-cs-faculty.stanford.edu/~uno/programs/commafree-eastman.w
  2. Recognize that the programming language used is C. Hence we need to use cweave and ctangle; not the original weave and tangle designed for Pascal.
    All those programs come bundled with TeX distributions like TeXLive. Once you have them installed and your system can find executables like pdflatex, it will also find cweave, etc.
  3. First, we generate the documentation.
    1. We run

      cweave commafree-eastman.w

      On my system this gives the output

      This is CWEAVE, Version 3.64 (TeX Live 2015/dev/Debian)
      *1*3*8
      Writing the output file...*1*3*8
      Writing the index...
      Done.
      (No errors were found.)
    2. Three files are generated, namely commafree-eastman.idx, commafree-eastman.scn and commafree-eastman.tex.
    3. As a small remark, the TeX file uses the cwebmac definitions. So it is a file full of definitions and this file is also provided with your TeX distribution:
      % head -n 1 commafree-eastman.tex
      \input cwebmac
      % locate cwebmac
      /usr/share/texlive/texmf-dist/tex/plain/cweb/cwebmac.tex
      /usr/share/texlive/texmf-dist/tex/plain/cweb/pdfXcwebmac.tex
      /usr/share/texlive/texmf-dist/tex/plain/cweb/pdfcwebmac.tex
      /usr/share/texlive/texmf-dist/tex/plain/cweb/pdfdcwebmac.tex
      /usr/share/texlive/texmf-dist/tex/plain/cweb/pdffcwebmac.tex
      /usr/share/texlive/texmf-dist/tex/plain/cweb/pdficwebmac.tex
    4. Then we run the original tex program to build a dvi file.

      tex commafree-eastman.tex

      This gives me the output:

      This is TeX, Version 3.14159265 (TeX Live 2015/dev/Debian) (preloaded format=tex)
      (./commafree-eastman.tex
      (/usr/share/texlive/texmf-dist/tex/plain/cweb/cwebmac.tex) *1 [1] *3 [2]
      [3] *8 Index: (./commafree-eastman.idx) [4] Section names:
      (./commafree-eastman.scn) [5] Table of contents: (./commafree-eastman.toc)
      [0] )
      Output written on commafree-eastman.dvi (6 pages, 16116 bytes).
      Transcript written on commafree-eastman.log.
    5. Followingly, a file commafree-eastman.dvi is written as the output states. Opening it directly with evince 3.14.2 shows the output with some errors (the typical font encoding problems of TeX like < in #include statements is replaced by ¡). With okular 0.21.3 I can view the document without any problems.
    6. We convert the DVI file using the utility dvipdfm provided with your TeX distribution:
      dvipdfm commafree-eastman.dvi

      For me, this gives:

      commafree-eastman.dvi -> commafree-eastman.pdf
      [1][2][3][4][5][6]
      45807 bytes written

      Then we have a PDF output file commafree-eastman.pdf containing the documentation.
      I think this is the usual way to view documents, you were looking for.

  4. We generate the source code:

    1. First we run ctangle:
      ctangle commafree-eastman.w

      which gives me:

      This is CTANGLE, Version 3.64 (TeX Live 2015/dev/Debian)
      *1*3*8
      Writing the output file (commafree-eastman.c):
      Done.
      (No errors were found.)
    2. Nice! So we got a commafree-eastman.c file.
    3. We compile the C file with a C compile like the GCC:
      gcc -Wall -o commafree-eastman commafree-eastman.c

      This gives me:

      ./commafree-eastman.w:38:1: warning: return type defaults to ‘int’ [-Wreturn-type]
       main (int argc,char*argv[]) {
       ^
      ./commafree-eastman.w: In function ‘main’:
      ./commafree-eastman.w:42:1: warning: control reaches end of non-void function [-Wreturn-type]
       }
       ^
    4. We run this new executable program:
      ./commafree-eastman
      Usage: ./commafree-eastman x1 x2 ... xn

      This gives the same output as in the video.

So we got a commafree-eastman.pdf and commafree-eastman.c file which was what we were going for.
I will also refer here to my Literate programming talk at PyGraz back then.

Running Donald Knuth’s CWEB programs

Mapping lexicographical comparison to integer comparison

Claim

(a, b) < (c, d) with a,b,c,d ∈ [0,M) ⇔ a · M + b < c · M + d
(a, b) = (c, d) with a,b,c,d ∈ [0,M) ⇔ a · M + b = c · M + d
(a, b) > (c, d) with a,b,c,d ∈ [0,M) ⇔ a · M + b > c · M + d

Proof by case distinction
relation to prove constraint evaluation conclusion
(a, b) < (c, d) a < c a · M + M – 1 < c · M
(a, b) < (c, d) a = c ∧ b < d a · M + b < c · M + d 0 + b < 0 + d ✓
(a, b) = (c, d) a = c ∧ b = d a · M + b = c · M + d 0 + 0 = 0 + 0 ✓
(a, b) > (c, d) a = c ∧ b > d a · M + b > c · M + d 0 + b > 0 + d ✓
(a, b) > (c, d) a > c a · M > c · M + M – 1
Rationale

This just exploits the nature of a lexicographical order itself for bounded domains. Consider two binary numbers and they get compared. This integer comparison equals the lexicographical tuple comparison of its digit expansion. I just wanted to make this explicit as my head was screwing around with it.

Mapping lexicographical comparison to integer comparison

LaTeX math cheatsheets

My girlfriend asked for LaTeX cheatsheets. There are many of them online, but she was specifically looking for notational commands for algebra and number theory. You need to be aware that notation is very specific for your teacher. So standardization has not been successful to raise mathematical symbols from a representational to a semantical level. In US/UK/AU (in my experience), there is little information encoded in notation and information is provided in short, concise sentences. In Europe, this is different and every teacher has his/her own notation.

But let’s get back to the point: I recommend the following cheatsheets, but evaluate yourself whether it fits your field.

If you need a more comprehensive reference, the LaTeX wikibook is a nice start. For a german audience, LaTeX@TUGraz also does a good job IMHO.

And as always: If you look for a symbol and have internet connection, detexify makes your life much easier!

LaTeX math cheatsheets

My second Kyusho experience

Today, on 15th of November 2015, our level 1 Kyusho course ended. I had my first experience with Kyusho in February when instructor Robert Göslbauer joined us in an introductory course for our association. Unlike the first time which was unstructured and rather application-oriented, this time he specifically covered all aspects of the Kyusho level 1 exam for all ~30 participants.

Robert showed us the majority of vital points along meridians described in Chinese medicine. They are used in a negative way potentially knocking out the opponent or losing neurological control in a certain part of the body. This demonstration uses vital point “Triple Warmer 17”, which can be very painful when applied properly. Vital points on the head like that one can make you feel dizzy, which is why regeneration can be important as can be seen in this video.

And why do I point out the word “experience” so much? I recognized two interesting effects: Even though you don’t move a lot, just by using the vital points, you feel exhausted afterwards. The pain itself affects your body similar to Budo training. You become more and more sensitive on the points. After applying the techniques multiple times, a little touch to the skin might already trigger pain (well, this is what my arm in the elbow felt like). However, this makes applying vital points in Jitsu combat more difficult, because in Jitsu (unlike Do arts like Aikido) you actually want to hurt the opponent and the first hit must be more precise than the following ones.

I am looking forward to taking the level 1 exam in the next months.

My second Kyusho experience

Cheatsheet for 1.0.0 releases

In the past 3 years I have released several software packages and a few of them exceeded the 1.0.0 release. I am not confident with it. I need structure. I need a plan for 1.0.0 releases.

So what is the background? In software engineering the 1.0.0 release is a famous one. The release 1.0.0 is considered as the first public release where people can depend on the API/ABI.

Consider a car. With the public release a car seat customizer will look at the measurements and publish his customized seats for this type of car. The car company should not publish the same car with different seat measurements a month later. This would destroy the business for the customizer and annoy costumers. This is about standardization.

The same applies to software engineering. With the 1.0.0 release you are expected not to change the components in such a way that depending software breaks. Like the car seat measurements. If you stick to semantic versioning as defined by semver.org, you are actually allowed to make “incompatible API changes” if you release by incrementing the first of three numbers. Like python was backwards incompatible when going from the 2.x to 3.x release. But semantic versioning is not yet common enough in the industry. So I stick to the paradigm that the 1.0.0 release should be well-considered.

And this is my personal check list for 1.0.0 releases in the future:

  1. I will define a target audience and consider the social impact.
  2. I have a documentation for the public API accessible in the web.
  3. I have a working testsuite.
  4. I have written several programs using the API and concluded that the API is convienient to use with the languages’ builtin.
  5. I will check my API against best practices in the software paradigm I am using.
  6. I have money & time to provide hotfixes 24 hours after release at the latest.
  7. I will specify
    1. metadata such as: contributors, license, version, date, contact options, associated URLs
    2. dependencies and system requirements
    3. simple installation instructions for the most common usecase
    4. performance estimates
    5. target audience, pitfalls & limitations
    6. a changelog (at least for major and minor releases)

Iff all these requirements are met, go ahead and release it as 1.0.0, but not otherwise.

Cheatsheet for 1.0.0 releases

Python riddle

Python 3.4.3 (default, Mar 26 2015, 22:03:40) 
[GCC 4.9.2] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> a = 42
>>> def f():
...     x = a
...     a = 3
... 
>>> f()
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "<stdin>", line 2, in f
UnboundLocalError: local variable 'a' referenced before assignment

Now let’s skip the second row of f() and do not define “a = 3”. Guess, what is the result:

>>> a = 42
>>> def f():
...     x = a
... 
>>> f()
>>>

This works the same in python 2 and 3. To me it is non-intuitive behavior, but gives reason why reusing variable names is generally a bad idea.

Python riddle