특정 속성을 가진 정수 찾기 - 프로젝트 오일러 문제 221
-
23-08-2019 - |
문제
나는 최근 프로젝트 Euler에 중독되어 노력하고 있습니다. 이것 다음! 나는 그것에 대한 분석을 시작했고 이미 문제를 크게 줄였습니다. 내 작업은 다음과 같습니다.
A = PQR 및
1/a = 1/p + 1/q + 1/r 그래서 pqr/a = pq + pr + qr
그리고 첫 번째 방정식 때문에 :
PQ+PR+QR = 1
정확히 p, q 및 r의 정확한 두 가지가 음수이어야하므로 다음 방정식을 찾는 것까지 단순화 할 수 있습니다.
ABC에 대한 AB = AC+BC+1
우리가 얻는 C를위한 해결 :
AB-1 = (A+B) c
C = (AB-1)/(A+B)
이것은 우리가 다음을 찾아야한다는 것을 의미합니다.
AB = 1 (모드 A+B)
그리고 A와 B와의 값은 다음과 같습니다.
a = abc = ab (ab-1)/(a+b)
많은 수학이라면 죄송합니다! 그러나 이제 우리가 다루어야 할 것은 하나의 조건과 두 방정식입니다. 이제 AB = 1 (모드 A+B)과 함께 AB (AB-1)/(A+B)로 작성된 150,000 번째 작은 정수를 찾아야하므로 이상적으로는 A IS를 검색하고 싶습니다. 가능한 한 작습니다.
쉽게 A <B를 가정하고 GCD (A, B) = 1을 알았습니다.
나의 첫 번째 구현은 간단하며 심지어 150,000 개의 솔루션을 충분히 빠르게 찾습니다. 그러나 150,000을 찾는 데 너무 오래 걸립니다. 가장 작은 솔루션. 어쨌든 코드는 다음과 같습니다.
n = 150000
seen = set()
a = 3
while len(seen) < n:
for b in range(2, a):
if (a*b)%(a+b) != 1: continue
seen.add(a*b*(a*b-1)//(a+b))
print(len(seen), (a, b), a*b*(a*b-1)//(a+b))
a += 1
다음 생각은 Stern-Brocot 나무를 사용하는 것이었지만 솔루션을 찾기에는 너무 느립니다. 최종 알고리즘은 중국의 나머지 정리를 사용하여 A+B의 다른 값이 수율 솔루션을 확인하는 것이 었습니다. 그 코드는 복잡하고 빠르지 만 빠르지는 않지만 충분하지 않습니다 ...
그래서 나는 절대적으로 아이디어가 없습니다! 누구든지 아이디어가 있습니까?
해결책
많은 프로젝트 Euler 문제와 마찬가지로, 비결은 무차별 인력 솔루션을보다 간단한 것으로 줄이는 기술을 찾는 것입니다.
A = pqr and
1/A = 1/p + 1/q + 1/r
그래서,
pq + qr + rp = 1 or -r = (pq - 1)/(p + q)
일반적인 손실없이 0 <p <-q <-r
k, 1 <= k <= p가 있습니다.
-q = p + k
-r = (-p(p + k) – 1) / (p + -p – k) = (p^2 + 1)/k + p
그러나 r은 정수이므로 k는 p^2 + 1을 나눕니다.
pqr = p(p + q)((p^2 + 1)/k + p)
따라서 A를 계산하려면 P를 반복해야하며 여기서 k는 P Squared Plus 1의 제수 인 값 만 가져갈 수 있습니다.
각 솔루션을 세트에 추가하면 필요한 150000th Alexandrian Integer를 찾을 때 중지 할 수 있습니다.
다른 팁
중국의 나머지, 빠른 구현에 관한이 기사는 도움이 될 수 있습니다 : www.codeproject.com/kb/recipes/crp.aspx
이것은 도구 및 라이브러리에 대한 더 많은 링크입니다.
도구 :
막시마http://maxima.sourceforge.net/Maxima는 차별화, 통합, Taylor 시리즈, Laplace 변환, 일반 미분 방정식, 선형 방정식, 다항식 및 세트, 목록, 벡터, 행렬 및 텐서를 포함한 기호 및 수치 표현 조작을위한 시스템입니다. Maxima는 정확한 분수, 임의의 정밀 정수 및 가변 정밀 부동 소수점 번호를 사용하여 높은 정밀 숫자 결과를 산출합니다. Maxima는 기능과 데이터를 2 차원에서 3 차원으로 플롯 할 수 있습니다.
수학http://mathomatic.org/math/Mathomatic은 무료, 휴대용 일반 목적 CAS (컴퓨터 대수 시스템) 및 방정식을 상징적으로 해결, 단순화, 결합 및 비교하고 복소수 및 다항식 산술 등을 수행 할 수있는 계산기 소프트웨어입니다. 사용.
scilab www.scilab.org/download/index_download.php scilab은 matlab 또는 simulink와 비슷한 수치 계산 시스템입니다. Scilab에는 수백 개의 수학적 기능이 포함되어 있으며 다양한 언어 (예 : C 또는 Fortran)의 프로그램이 대화식으로 추가 될 수 있습니다.
Mathstudio Mathstudio.sourceforge.net 대화식 방정식 편집기 및 단계별 솔버.
도서관:
Armadillo C ++ 라이브러리http://arma.sourceforge.net/Armadillo C ++ 라이브러리는 간단하고 사용하기 쉬운 인터페이스를 갖는 동시에 선형 대수 작업 (매트릭스 및 벡터 수학)을위한 효율적인 기반을 제공하는 것을 목표로합니다.
Blitz ++http://www.oonumerics.org/blitz/ Blitz ++는 과학 컴퓨팅을위한 C ++ 클래스 라이브러리입니다.
biginteger c#http://msdn.microsoft.com/pt-br/magazine/cc163441.aspx
libapmathhttp://freshmeat.net/projects/libapmathApmath-Project의 홈페이지에 오신 것을 환영합니다. 이 프로젝트의 목표는 임의의 정밀 C ++-라이브러리를 구현하는 것입니다.이 라이브러리는 가장 편리합니다. 이는 모든 작업이 운영자-과부하로 구현되며, 이름 지정은 대부분의 것과 동일하다는 것을 의미합니다.
libmathttp://freshmeat.net/projects/libmatMAT는 C ++ 수학 템플릿 클래스 라이브러리입니다. 다양한 매트릭스 작업 에이 라이브러리를 사용하고 다항식의 루트 찾기, 방정식 해결 등을 사용하십시오. 라이브러리에는 C ++ 헤더 파일 만 포함되므로 컴파일이 필요하지 않습니다.
애니 마스http://www.yonsen.bz/animath/animath.htmlAnimath는 C ++로 완전히 구현 된 유한 요소 방법 라이브러리입니다. 유체 구조 상호 작용 시뮬레이션에 적합하며 수학적으로 고차 사면체 요소를 기반으로합니다.
확인. 다음은 내 중국의 나머지 정리 솔루션으로 더 많은 연주를합니다. A+B P = 1 (모드 4). 이것은 2, 5, 13, 17, 29, 37과 같은 프라임의 배수 인 A+B 만 확인하면 더 빠른 계산을 허용합니다.
따라서 가능한 A+B 값의 샘플이 있습니다.
[5, 8, 10, 13, 16, 17, 20, 25, 26, 29, 32, 34, 37, 40, 41, 50, 52, 53, 58, 61, 64, 65, 68, 73, 74, 80, 82, 85, 89, 97, 100]
그리고 여기에 중국의 나머지 정리를 사용한 전체 프로그램이 있습니다.
cachef = {}
def factors(n):
if n in cachef: cachef[n]
i = 2
while i*i <= n:
if n%i == 0:
r = set([i])|factors(n//i)
cachef[n] = r
return r
i += 1
r = set([n])
cachef[n] = r
return r
cachet = {}
def table(n):
if n == 2: return 1
if n%4 != 1: return
if n in cachet: return cachet[n]
a1 = n-1
for a in range(1, n//2+1):
if (a*a)%n == a1:
cachet[n] = a
return a
cacheg = {}
def extended(a, b):
if a%b == 0:
return (0, 1)
else:
if (a, b) in cacheg: return cacheg[(a, b)]
x, y = extended(b, a%b)
x, y = y, x-y*(a//b)
cacheg[(a, b)] = (x, y)
return (x, y)
def go(n):
f = [a for a in factors(n)]
m = [table(a) for a in f]
N = 1
for a in f: N *= a
x = 0
for i in range(len(f)):
if not m[i]: return 0
s, t = extended(f[i], N//f[i])
x += t*m[i]*N//f[i]
x %= N
a = x
while a < n:
b = n-a
if (a*b-1)%(a+b) == 0: return a*b*(a*b-1)//(a+b)
a += N
li = [5, 8, 10, 13, 16, 17, 20, 25, 26, 29, 32, 34, 37, 40, 41, 50, 52, 53, 58, 61, 64, 65, 68, 73, 74, 80, 82, 85, 89, 97, 100]
r = set([6])
find = 6
for a in li:
g = go(a)
if g:
r.add(g)
#print(g)
else:
pass#print(a)
r = list(r)
r.sort()
print(r)
print(len(r), 'with', len(li), 'iterations')
이것은 더 좋지만 더 개선하기를 희망합니다 (예 : A+B = 2^N은 결코 해결책이 아닌 것 같습니다).
또한 다음과 같은 기본 대체를 고려하기 시작했습니다.
a = u+1 및 b = v+1
AB = 1 (모드 A+B)
UV+U+V = 0 (모드 U+V+2)
그러나 나는 그것에 대해 많은 개선을 볼 수 없습니다 ...