一点一点前进...

0%

欧拉项目 | 518题 | 质数三元组和等比数列

发表于 分类于 ， Disqus：

518题链接

S(n) = sum(a+b+c)，其中a b c满足以下三个性质：

  1. a b c都是质数
  2. a < b < c < n
  3. a+1, b+1, c+1能组成等比数列

题目中给出了一个例子S(100) = 1035可以用于测试。
求S(10^8)

前一亿个数里面有576万个质数，显然，不可以接受O(n^2)的算法。
为了描述方便，后面a b c是题目中的a+1 b+1和c+1。
a b c能组成等比数列，那么b^2=a*c。
给定一个a，我们不能用遍历的方式去找b，否则会是O(n^2)的算法。我们需要一个更快的方式找到可能的b。
考虑这样一个问题，a的质因数里面有4个p，那么b的质因数里面至少有2个p，那么a和b之间最少也差着p^2；a的质因数里面有5个q，那么b的质因数里面至少有3个q，那么a和b之间最少也差着q^3；同时考虑p和q的话，那么差值是p^2*q^3。有了这个差值，就能够大大减少内层循环的数目。通过计数器发现，内层循环大概是8-9次的。
问题转换到如何得到a的质因数及其次数，这是一个比较简单的问题。 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
public static List<Tuple<long, int>> TrialDivisioFactor(long n, long[] primes)
{
    var results = new List<Tuple<long, int>>();
    int index = 0;
    while (true)
    {
        if (primes[index] != 0)
        {
            if (primes[n] != 0)
            {
                results.Add(new Tuple<long, int>(n, 1));
                break;
            }

            int count = 0;
            while (n % primes[index] == 0)
            {
                count++;
                n /= primes[index];
            }

            if (count != 0)
            {
                results.Add(new Tuple<long, int>(primes[index], count));
            }

            if (n == 1)
            {
                break;
            }
        }

        index++;
    }

    return results;
}

GenPrimes就是通过筛选法得到质数。
下面是完整的程序。 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
int end = 100000000;
long sum = 0;
Utils.GenPrimes(end, out long[] primes);
for (int i = 0; i < primes.Length; i++)
{
    if (primes[i] != 0)
    {
        var a = primes[i] + 1;
        var factors = Utils.TrialDivisioFactor(a, primes);
        var step = 1;
        foreach (var f in factors)
        {
            step *= (int)Math.Pow(f.Item1, (f.Item2 + 1) / 2);
        }
        for (long b = a + step; b < end; b += step)
        {
            var c = b * b / a;
            if (c >= end)
            {
                break;
            }

            if (primes[b - 1] != 0 && primes[c - 1] != 0)
            {
                sum += a + b + c - 3;
            }
        }
    }
}

return sum;