Find the largest gap between good primes Following the fine tradition of questions such as Find the largest prime whose length, sum and product is prime , this is a variant on a largest prime challenge.
Input
Your code should not take any input.
Definition
We say a prime p is good if p-1 has exactly 2 distinct prime factors.
Output
Your code should output the absolute difference between consecutive good primes q and p so that |q-p| is as large as possible and q is the smallest good prime larger than p.  You may output any number of good pairs and your last output will be taken as the score.
Example
The sequence of the first 55 good primes  is https://oeis.org/A067466 .
Score
Your score is simply |q-p| for the pair of good primes you output.
Languages and libraries
You can use any language or library you like (that wasn't designed for this challenge) except for any library functions for primality testing or factoring integers. However, for the purposes of scoring I will run your code on my machine so please provide clear instructions for how to run it on Ubuntu.
My Machine The timings will be run on my machine. This is a standard Ubuntu install on an 8GB AMD FX-8350 Eight-Core Processor.  This also means I need to be able to run your code.
Details

I will kill your code after 2 minutes unless it starts to run out of memory before that. It should therefore make sure to output something before the cut off.
You may not use any external source of primes.
You may use probabilistic prime testing methods although I am told by Mego that with good tables, Miller-Rabin can test up to 341,550,071,728,321 (or even higher) deterministically. See also http://miller-rabin.appspot.com/ .

Best entries that check all integers from 1

756 by cat in Go
756 by El'endia Starman in Python
1932 by Adnan in C# (using mono 3.2.8)
2640 by yeti in Python (using pypy 4.01)
2754 by Reto Koradi in C++
3486 by Peter Taylor in Java
3900 by primo in RPython (using pypy 4.01)
4176 by The Coder in Java

Best entries that may skip a large number of integers to find a large gap

14226 by Reto Koradi in C++
22596 by primo in RPython (using pypy 4.01). Record reached after 5 seconds!
Probably 4032, mixed Atkin-Bernstein sieve and "deterministic" Miller-Rabin

Wheel factorisation and good primes

It's highly obvious that with the exceptions of 2, 3, and 5, every prime is coprime to 2*3*5 = 60. There are 16 equivalence classes modulo 60 which are coprime to 60, so any primality test only needs to check those 16 cases.

However, when we're looking for "good" primes we can thin the herd even more. If gcd(x, 60) = 1, we find that in most cases gcd(x-1, 60) is either 6 or 10. There are 6 values of x for which it is 2:

17, 23, 29, 47, 53, 59


Therefore we can precompute the "good" primes of the form 2^a 3^b + 1 and 2^a 5^b + 1 and merge them into the result of a prime generator which only considers 10% of numbers as even potential primes.

Implementation notes

Since I already had a Java implementation of the Atkin-Bernstein sieve lying around, and that already uses a wheel as a key component, it seemed natural to strip out the unnecessary spokes and adapt the code. I originally tried using a producer-consumer architecture to exploit the 8 cores, but memory management was too messy.

To test whether a prime is a "good" prime, I'm using a "deterministic" Miller-Rabin test (which really means a Miller-Rabin test which someone else has pre-checked against a list generated deterministically). This can certainly be rewritten to also use Atkin-Bernstein, with some caching to cover the ranges corresponding to sqrt, cbrt, etc., but I'm not sure whether it would be an improvement (because it would be testing lots of numbers which I don't need to test), so that's an experiment for another day.

On my fairly old computer this runs to 

987417437 - 987413849 = 3588
1123404923 - 1123401023 = 3900
1196634239 - 1196630297 = 3942
1247118179 - 1247114147 = 4032


in pretty much exactly 2 minutes, and then

1964330609 - 1964326433 = 4176
2055062753 - 2055058529 = 4224
2160258917 - 2160254627 = 4290


at 3:10, 3:20, and 3:30 respectively.

import java.util.*;

public class PPCG65876 {
    public static void main(String[] args) {
        long[] specials = genSpecials();
        int nextSpecialIdx = 0;
        long nextSpecial = specials[nextSpecialIdx];
        long p = 59;
        long bestGap = 2;

        for (long L = 1; true; L += B) {

            long[][] buf = new long[6][B &gt;&gt; 6];
            int[] Lmodqq = new int[qqtab.length];
            for (int i = 0; i &lt; Lmodqq.length; i++) Lmodqq[i] = (int)(L % qqtab[i]);

            for (long[] arr : buf) Arrays.fill(arr, -1); // TODO Can probably get a minor optimisation by inverting this
            for (int[] parms : elliptic) traceElliptic(buf[parms[0]], parms[1], parms[2], parms[3] - L, parms[4], parms[5], Lmodqq, totients[parms[0]]);
            for (int[] parms : hyperbolic) traceHyperbolic(buf[parms[0]], parms[1], parms[2], parms[3] - L, Lmodqq, totients[parms[0]]);

            // We need to filter down to squarefree numbers.
            long pg_base = L * M;
            squarefreeMid(buf, invTotients, pg_base, 247, 1, 38);
            squarefreeMid(buf, invTotients, pg_base, 253, 1, 38);
            squarefreeMid(buf, invTotients, pg_base, 257, 1, 38);
            squarefreeMid(buf, invTotients, pg_base, 263, 1, 38);
            squarefreeMid(buf, invTotients, pg_base, 241, 0, 2);
            squarefreeMid(buf, invTotients, pg_base, 251, 0, 2);
            squarefreeMid(buf, invTotients, pg_base, 259, 0, 2);
            squarefreeMid(buf, invTotients, pg_base, 269, 0, 2);

            // Turn bitmasks into primes
            long[] page = new long[150000]; // TODO This can almost certainly be smaller
            long[] transpose = new long[6];
            for (int j = 0, off = 0; j &lt; (B &gt;&gt; 6); j++) {
                // Reduce cache locality problems by transposing.
                for (int k = 0; k &lt; 6; k++) transpose[k] = buf[k][j];
                for (long mask = 1L; mask != 0; mask &lt;&lt;= 1) {
                    for (int k = 0; k &lt; 6; k++) {
                        if ((transpose[k] &amp; mask) == 0) page[off++] = pg_base + totients[k];
                    }

                    pg_base += M;
                }
            }

            // Insert specials and look for gaps.
            for (long q : page) {
                if (q == 0) break;

                // Do we need to insert one or more specials?
                while (nextSpecial &lt; q) {
                    if (nextSpecial - p &gt; bestGap) {
                        bestGap = nextSpecial - p;
                        System.out.format("%d - %d = %d\n", nextSpecial, p, bestGap);
                    }

                    p = nextSpecial;
                    nextSpecial = specials[++nextSpecialIdx];
                }

                if (isGood(q)) {
                    if (q - p &gt; bestGap) {
                        bestGap = q - p;
                        System.out.format("%d - %d = %d\n", q, p, bestGap);
                    }

                    p = q;
                }
            }

        }
    }

    static long[] genSpecials() {
        // 2^a 3^b + 1 or 2^a 5^b + 1
        List&lt;Long&gt; tmp = new LinkedList&lt;Long&gt;();
        for (long threes = 3; threes &lt;= 4052555153018976267L; threes *= 3) {
            for (long t = threes; t &gt; 0; t &lt;&lt;= 1) tmp.add(t + 1);
        }
        for (long fives = 5; fives &lt;= 7450580596923828125L; fives *= 5) {
            for (long f = fives; f &gt; 0; f &lt;&lt;= 1) tmp.add(f + 1);
        }

        // Filter down to primes
        Iterator&lt;Long&gt; it = tmp.iterator();
        while (it.hasNext()) {
            long next = it.next();
            if (next &lt; 60 || next &gt; 341550071728321L || !isPrime(next)) it.remove();
        }

        Collections.sort(tmp);
        long[] specials = new long[tmp.size()];
        for (int i = 0; i &lt; tmp.size(); i++) specials[i] = tmp.get(i);

        return specials;
    }

    private static boolean isGood(long p) {
        long d = p - 1;
        while ((d &amp; 1) == 0) d &gt;&gt;= 1;

        if (d == 1) return false;

        // Is d a prime power?
        if (d % 3 == 0 || d % 5 == 0) {
            // Because of the way the filters before this one work, nothing should reach here.
            throw new RuntimeException("Should be unreachable");
        }

        // TODO Is it preferable to reuse the Atkin-Bernstein code, caching pages which correspond
        // to the possible power candidates?
        if (isPrime(d)) return true;
        for (int a = (d % 60 == 1 || d % 60 == 49) ? 2 : 3; (1L &lt;&lt; a) &lt; d; a++) {
            long r = (long)(0.5 + Math.pow(d, 1. / a));
            if (d == (long)(0.5 + Math.pow(r, a)) &amp;&amp; isPrime(r)) return true;
        }

        return false;
    }

    /*---------------------------------------------------
               Deterministic Miller-Rabin
    ---------------------------------------------------*/
    public static boolean isPrime(int x) {
        // See isPrime(long). We pick bases which are known to work for the entire range of int.
        // Special case for the bases.
        if (x == 2 || x == 7 || x == 61) return true;

        int d = x - 1;
        int s = 0;
        while ((d &amp; 1) == 0) { s++; d &gt;&gt;= 1; } // TODO Can be optimised

        if (!isSPRP(2, d, s, x)) return false;
        if (!isSPRP(7, d, s, x)) return false;
        if (!isSPRP(61, d, s, x)) return false;
        return true;
    }

    private static boolean isSPRP(int b, int d, int s, int x /* == d &lt;&lt; s */) {
        int l = modPow(b, d, x);
        if (l == 1 || l == x - 1) return true;
        for (int r = 1; r &lt; s; r++) {
            l = modPow(l, 2, x);
            if (l == x - 1) return true;
            if (l == 1) return false;
        }

        return false;
    }

    public static int modPow(int a, int b, int c) {
        int accum = 1;
        while (b &gt; 0) {
            if ((b &amp; 1) == 1) accum = (int)(accum * (long)a % c);
            a = (int)(a * (long)a % c);
            b &gt;&gt;= 1;
        }
        return accum;
    }

    public static boolean isPrime(long x) {
        if (x &lt; Integer.MAX_VALUE) return isPrime((int)x);

        long d = x - 1;
        int s = 0;
        while ((d &amp; 1) == 0) { s++; d &gt;&gt;= 1; } // TODO Can be optimised

        // If b^d == 1 (mod x) or (b^d)^(2^r) == -1 (mod x) for some r &lt; s then we pass for base b.
        // We select bases according to Jaeschke, Gerhard (1993), "On strong pseudoprimes to several bases", Mathematics of Computation 61 (204): 915–926, doi:10.2307/2153262
        // TODO Would it be better to use a set of 5 bases from http://miller-rabin.appspot.com/ ?
        if (!isSPRP(2, d, s, x)) return false;
        if (!isSPRP(3, d, s, x)) return false;
        if (!isSPRP(5, d, s, x)) return false;
        if (!isSPRP(7, d, s, x)) return false;
        if (x &lt; 3215031751L) return true;
        if (!isSPRP(11, d, s, x)) return false;
        if (x &lt; 2152302898747L) return true;
        if (!isSPRP(13, d, s, x)) return false;
        if (x &lt; 3474749660383L) return true;
        if (!isSPRP(17, d, s, x)) return false;
        if (x &lt; 341550071728321L) return true;

        throw new IllegalArgumentException("Overflow");
    }

    private static boolean isSPRP(long b, long d, int s, long x /* == d &lt;&lt; s */) {
        if (b * (double)x &gt; Long.MAX_VALUE) throw new IllegalArgumentException("Overflow"); // TODO Work out more precise page bounds

        long l = modPow(b, d, x);
        if (l == 1 || l == x - 1) return true;
        for (int r = 1; r &lt; s; r++) {
            l = modPow(l, 2, x);
            if (l == x - 1) return true;
            if (l == 1) return false;
        }

        return false;
    }

    /**
     * Computes a^b (mod c). We assume c &amp;lt; 2^62.
     */
    public static long modPow(long a, long b, long c) {
        long accum = 1;
        while (b &gt; 0) {
            if ((b &amp; 1) == 1) accum = prodMod(accum, a, c);
            a = prodMod(a, a, c);
            b &gt;&gt;= 1;
        }
        return accum;
    }

    /**
     * Computes a*b (mod c). We assume c &amp;lt; 2^62.
     */
    private static long prodMod(long a, long b, long c) {
        // The naive product would require 128-bit integers.

        // Consider a = (A &lt;&lt; 32) + B, b = (C &lt;&lt; 31) + D. Different shifts chosen deliberately.
        // Then ab = (AC &lt;&lt; 63) + (AD &lt;&lt; 32) + (BC &lt;&lt; 31) + BD with intermediate values remaining in 63 bits.
        long AC = (a &gt;&gt; 32) * (b &gt;&gt; 31) % c;
        long AD = (a &gt;&gt; 32) * (b &amp; ((1L &lt;&lt; 31) - 1)) % c;
        long BC = (a &amp; ((1L &lt;&lt; 32) - 1)) * (b &gt;&gt; 31) % c;
        long BD = (a &amp; ((1L &lt;&lt; 32) - 1)) * (b &amp; ((1L &lt;&lt; 31) - 1)) % c;

        long t = AC;
        for (int i = 0; i &lt; 31; i++) {
            t = (t + t) % c;
        }
        // t = (AC &lt;&lt; 31)
        t = (t + AD) % c;
        t = (t + t) % c;
        t = (t + BC) % c;
        // t = (AC &lt;&lt; 32) + (AD &lt;&lt; 1) + BC
        for (int i = 0; i &lt; 31; i++) {
            t = (t + t) % c;
        }
        // t = (AC &lt;&lt; 63) + (AD &lt;&lt; 32) + (BC &lt;&lt; 31)
        return (t + BD) % c;
    }

    /*---------------------------------------------------
                      Atkin-Bernstein
    ---------------------------------------------------*/
    // Page size.
    private static final int B = 1001 &lt;&lt; 6;
    // Wheel modulus for sharding between binary quadratic forms.
    private static final int M = 60;

    // Squares of primes 5 &lt; q &lt; 240
    private static final int[] qqtab = new int[] {
        49, 121, 169, 289, 361, 529, 841, 961, 1369, 1681, 1849, 2209, 2809,
        3481, 3721, 4489, 5041, 5329, 6241, 6889, 7921, 9409, 10201, 10609, 11449, 11881,
        12769, 16129, 17161, 18769, 19321, 22201, 22801, 24649, 26569, 27889, 29929, 32041, 32761,
        36481, 37249, 38809, 39601, 44521, 49729, 51529, 52441, 54289, 57121
    };
    // If a_i == q^{-2} (mod 60) is the reciprocal of qq[i], qq60tab[i] = qq[i] + (1 - a_i * qq[i]) / 60
    private static int[] qq60tab = new int[] {
        9, 119, 31, 53, 355, 97, 827, 945, 251, 1653, 339, 405, 515,
        3423, 3659, 823, 4957, 977, 6137, 1263, 7789, 1725, 10031, 1945, 2099, 11683,
        2341, 2957, 16875, 3441, 18999, 21831, 22421, 4519, 4871, 5113, 5487, 31507, 32215,
        35873, 6829, 7115, 38941, 43779, 9117, 9447, 51567, 9953, 56169
    };

    /**
     * Produces a set of parameters for traceElliptic to find solutions to ax^2 + cy^2 == d (mod M).
     * @param d The target residue.
     * @param a Binary quadratic form parameter.
     * @param c Binary quadratic form parameter.
     */
    private static List&lt;int[]&gt; initElliptic(final int[] invTotients, final int d, final int a, final int c) {
        List&lt;int[]&gt; rv = new ArrayList&lt;int[]&gt;();

        // The basic idea is that we maintain an invariant of the form
        //     M k = a x^2 + c y^2 - d
        // Therefore we increment x in steps F such that
        //     a((x + F)^2 - x^2) == 0 (mod M)
        // and similarly for y in steps G.
        int F = computeIncrement(a, M), G = computeIncrement(c, M);
        for (int f = 1; f &lt;= F; f++) {
            for (int g = 1; g &lt;= G; g++) {
                if ((a*f*f + c*g*g - d) % M == 0) {
                    rv.add(new int[] { invTotients[d], (2*f + F)*a*F/M, (2*g + G)*c*G/M, (a*f*f + c*g*g - d)/M, 2*a*F*F/M, 2*c*G*G/M });
                }
            }
        }

        return rv;
    }

    private static int computeIncrement(int a, int M) {
        // Find smallest F such that M | 2aF and M | aF^2
        int l = M / gcd(M, 2 * a);
        for (int F = l; true; F += l) {
            if (a*F*F % M == 0) return F;
        }
    }

    public static int gcd(int a, int b) {
        while (b != 0) {
            int t = b;
            b = a % b;
            a = t;
        }

        return a;
    }

    // NB This is generalised somewhat from primegen's implementation.
    private static void traceElliptic(final long[] buf, int x, int y, long start, final int cF2, final int cG2, final int[] Lmodqq, final int d) {
        // Bring the annular segment into the range of ints.
        start += 1000000000;
        while (start &lt; 0) {
            start += x;
            x += cF2;
        }
        start -= 1000000000;
        int i = (int)start;

        while (i &lt; B) {
            i += x;
            x += cF2;
        }

        while (true) {
            x -= cF2;
            if (x &lt;= cF2 &gt;&gt; 1) {
                // It makes no sense that doing this in here should perform well, but empirically it does much better than
                // only eliminating the squares once.
                squarefreeTiny(buf, Lmodqq, d);
                return;
            }
            i -= x;

            while (i &lt; 0) {
                i += y;
                y += cG2;
            }

            int i0 = i, y0 = y;
            while (i &lt; B) {
                buf[i &gt;&gt; 6] ^= 1L &lt;&lt; i;
                i += y;
                y += cG2;
            }
            i = i0;
            y = y0;
        }
    }

    // This only handles 3x^2 - y^2, and is closer to a direct port of primegen.
    private static void traceHyperbolic(final long[] a, int x, int y, long start, final int[] Lmodqq, final int d) {
        x += 5;
        y += 15;

        // Bring the segment into the range of ints.
        start += 1000000000;
        while (start &lt; 0) {
            start += x;
            x += 10;
        }
        start -= 1000000000;
        int i = (int)start;

        while (i &lt; 0) {
            i += x;
            x += 10;
        }

        while (true) {
            x += 10;
            while (i &gt;= B) {
                if (x &lt;= y) {
                    squarefreeTiny(a, Lmodqq, d);
                    return;
                }
                i -= y;
                y += 30;
            }

            int i0 = i, y0 = y;
            while (i &gt;= 0 &amp;&amp; y &lt; x) {
                a[i &gt;&gt; 6] ^= 1L &lt;&lt; i;
                i -= y;
                y += 30;
            }
            i = i0 + x - 10;
            y = y0;
        }
    }

    private static void squarefreeTiny(final long[] a, final int[] Lmodqq, final int d) {
        for (int j = 0; j &lt; qqtab.length; ++j) {
            int qq = qqtab[j];
            int k = qq - 1 - ((Lmodqq[j] + qq60tab[j] * d - 1) % qq);
            while (k &lt; B) {
                a[k &gt;&gt; 6] |= 1L &lt;&lt; k;
                k += qq;
            }
        }
    }

    private static void squarefreeMid(long[][] buf, int[] invTotients, final long base, int q, int dqq, int di) {
        int qq = q * q;
        q = M * q + (M * M / 4);

        while (qq &lt; M * B) {
            int i = qq - (int)(base % qq);
            if ((i &amp; 1) == 0) i += qq;

            if (i &lt; M * B) {
                int qqhigh = ((qq / M) &lt;&lt; 1) + dqq;
                int ilow = i % M;
                int ihigh = i / M;
                while (ihigh &lt; B) {
                    int n = invTotients[ilow];
                    if (n &gt;= 0) buf[n][ihigh &gt;&gt; 6] |= 1L &lt;&lt; ihigh;

                    ilow += di;
                    ihigh += qqhigh;
                    if (ilow &gt;= M) {
                        ilow -= M;
                        ihigh += 1;
                    }
                }
            }

            qq += q;
            q += M * M / 2;
        }

        squarefreebig(buf, invTotients, base, q, qq);
    }

    private static void squarefreebig(long[][] buf, int[] invTotients, final long base, int q, long qq) {
        long bound = base + M * B;
        while (qq &lt; bound) {
            long i = qq - (base % qq);
            if ((i &amp; 1) == 0) i += qq;

            if (i &lt; M * B) {
                int pos = (int)i;
                int n = invTotients[pos % M];
                if (n &gt;= 0) {
                    int ihigh = pos / M;
                    buf[n][ihigh &gt;&gt; 6] |= 1L &lt;&lt; ihigh;
                }
            }

            qq += q;
            q += M * M / 2;
        }
    }

    // The relevant totients of M - those which only have one forced prime factor.
    static final int[] totients = new int[] { 17, 23, 29, 47, 53, 59 };
    private static final int[] invTotients;
    // Parameters for tracing the hyperbolic BQF used for 59+60Z.
    private static final int[][] hyperbolic = new int[][] {
        {5, 1, 2, -1}, {5, 1, 8, -2}, {5, 1, 22, -9}, {5, 1, 28, -14}, {5, 4, 7, -1}, {5, 4, 13, -3}, {5, 4, 17, -5}, {5, 4, 23, -9},
        {5, 5, 4, 0}, {5, 5, 14, -3}, {5, 5, 16, -4}, {5, 5, 26, -11}, {5, 6, 7, 0}, {5, 6, 13, -2}, {5, 6, 17, -4}, {5, 6, 23, -8},
        {5, 9, 2, 3}, {5, 9, 8, 2}, {5, 9, 22, -5}, {5, 9, 28, -10}, {5, 10, 1, 4}, {5, 10, 11, 2}, {5, 10, 19, -2}, {5, 10, 29, -10}
    };

    // Parameters for tracing the elliptic BQFs used for all totients except 11 and 59.
    private static final int[][] elliptic;
    static {
        invTotients = new int[M];
        Arrays.fill(invTotients, -1);
        for (int i = 0; i &lt; totients.length; i++) invTotients[totients[i]] = i;

        // Calculate the parameters for tracing the elliptic BQFs from a table of the BQF used for each totient.
        // E.g. for 17+60Z we use 5x^2 + 3y^2.
        int[][] bqfs = new int[][] {
            {17, 5, 3}, {23, 5, 3}, {29, 4, 1}, {47, 5, 3}, {53, 5, 3}
        };
        List&lt;int[]&gt; parmSets = new ArrayList&lt;int[]&gt;();
        for (int[] bqf : bqfs) parmSets.addAll(initElliptic(invTotients, bqf[0], bqf[1], bqf[2]));
        elliptic = parmSets.toArray(new int[0][]);
    }
}


Save as PPCG65876.java, compile as javac PPCG65876.java, and run as java -Xmx1G PPCG65876.