2

对于周长为 p 的整数直角三角形有许多解(a,b,c),对于所有这些解,a+b+c == p 和勾股定理也适用。我正在编写一个 Python 脚本来计算周长 <= 1000 的三角形可能的最大解决方案数。

我的脚本是正确的,但是它需要永远运行。我敢肯定即使使用我的 i7 处理器也需要 30 多分钟,所以我需要对其进行优化。有人能帮我吗?(这是 Project Euler 的一个问题,以防有人想知道)

def solutions(p):
    result = []

    for a in range(1, p + 1):
        for b in range(1, p - a + 1):
            for c in range(1, p - a - b + 1):
                if a + b + c == p and a < b and b < c:
                    d = a ** 2
                    e = b ** 2
                    f = c ** 2

                    if (d + e == f) or (e + f == d) or (f + d == e):
                        result.append((a, b, c))
    return len(result)


max_p = 0
max_solutions = 0

for p in range(3, 1001):
    print("Processing %d" % p)
    s = solutions(p)

    if s > max_solutions:
        max_solutions = s
        max_p = p

print("%d has %d solutions" % (max_p, max_solutions))
4

4 回答 4

1

一个更好的:

def solution(n):
    count = 0
    for c in range(n // 3 + 1, n // 2):
        for a in range(1, n // 3):
            b = n - a - c
            if b <= 0:
                continue
            if a >= b:
                continue
            if a * a + b * b != c * c:
                continue
            count += 1
    return count
于 2013-02-07T05:22:20.597 回答
1

这是我对您的程序的重写。

首先,我预先计算所有平方值。这不仅避免了乘法,而且意味着 Python 不会不断地为所有平方值创建和垃圾收集对象。

接下来,我摆脱了三角形第三边的循环。一旦您选择了值a并且b只有一个可能的值符合标准a + b + c == 1000,所以这只是测试那个。这将问题从大约 O(n^3) 变为大约 O(n^2),这是一个巨大的改进。

然后我尝试运行它。在我用了 4 年的电脑上,它在大约 46 秒内完成,所以我停在那里,你去吧。

我做了一个谷歌搜索,发现了这个问题的讨论;如果我看到的讨论是正确的,那么这个程序会打印正确的答案。

upper_bound = 1000

sqr = [i**2 for i in range(upper_bound+1)]

def solutions(p):
    result = []

    for a in range(1, p - 1):
        for b in range(1, p - a):
            c = p - (a + b)
            if a < b < c:
                d = sqr[a]
                e = sqr[b]
                f = sqr[c]

                if (d + e == f) or (e + f == d) or (f + d == e):
                    result.append((a, b, c))
    return len(result)


max_p = 0
max_solutions = 0

for p in range(3, upper_bound+1):
    print("Processing %d" % p)
    s = solutions(p)

    if s > max_solutions:
        max_solutions = s
        max_p = p

print("%d has %d solutions" % (max_p, max_solutions))

编辑:这是我正在玩的一个更快的版本。它在评论中包含了来自@gnibbler 的建议。

upper_bound = 1000

sqr = [i**2 for i in range(upper_bound+1)]

def solution(p):
    count = 0
    for a in range(1, p - 1):
        for b in range(a, p - a):
            c = p - (a + b)
            d = sqr[a]
            e = sqr[b]
            f = sqr[c]

            if (d + e == f):
                count += 1
    return count

c, p = max((solution(p), p) for p in range(3, upper_bound+1))
print("%d has %d solutions" % (p, c))

在我的电脑上,这个版本需要 31 秒而不是 46 秒。

使用的棘手业务max()并没有真正使它明显更快。我在没有预先计算平方的情况下尝试了它,它非常慢,太慢了,我不想等待确切的时间。

于 2013-02-07T05:24:07.467 回答
0

知道了。它只是依赖于设置 a^2+b^2=c^2 然后替换 p - a - b = c

 1 from math import pow
  2 
  3 def see_if_right_triangle(p):
        solutions = 0
  4     # Accepts the perimeter as input
  5     for a in range(1, p): 
  6         for b in range(1, p):
  7             if 2*p*b + 2*p*a - pow(p, 2) == 2*a*b:
  8                solutions += 1
        print "The perimeter {p} has {sol} number of solutions".format(p=p, sol=solutions)
 10                        
 11 
 12 for p in range(3, 1001):
 13     see_if_right_triangle(p)

我认为这可以进一步优化......特别是如果你找出一些数学来缩小你将接受的范围 a 和 b

于 2013-02-07T05:29:30.687 回答
0

这不是您优化的代码,而是我自己的代码(我用于此问题)。我开始做一些代数,以使程序非常简单,并且不必迭代1000^3次数(1-1000 用于a,然后 1-1000 用于的每个ba,1-1000 用于.cb

# Project Euler 9
'''
Algebra behind Final Method:
a + b + c = 1000 | a^2 + b^2 = c^2
c = 1000 - (a + b) # Solving for C
a^2 + b^2 = (1000 - (a + b))^2 # Substituting value for C
a^2 + b^2 = 1000000 - 2000a - 2000b + (a + b)^2 # simplifying
a^2 + b^2 = 1000000 - 2000a - 2000b + a^2 + b^2 + 2ab # simplifying again
0 = 1000000 - 2000a - 2000b + 2ab # new formula
2000a - 2ab = 1000000 - 2000b # isolating A
1000a - ab = 500000 - 1000b # divide by 2 to simplify
a(1000 - b) = 500000 - 1000b # factor out A
a = (500000 - 1000b) / (1000 - b) # solve for A
'''
def pE9():
   from math import sqrt
   a, b, c = 1, 1, 1
   while True:
      b += 1
      a = (500000 - 1000 * b) / (1000 - b)
      c = sqrt(a**2 + b**2)
      if a + b + c == 1000 and a**2 + b**2 == c**2:
         break
   print int(a * b * c)

from timeit import timeit
print timeit(pE9, number = 1)

使用number = 1所以只需测试一次。

输出:

>>> 
# Answer hidden
0.0142664994414
于 2013-02-07T18:49:48.307 回答