ディオファントス方程式
(1) 1/a+1/b=p/10 (2) 1/a+1/b=p/100 (3) 1/a+1/b=p/1000 (a,b,pは正の整数で、a≦b)
について、各方程式で、(a,b,p)の解が存在できない最小のpを、各場合について求めて下さい。
らすかるさんからのコメントです。(令和7年5月31日付け)
(1) 式を変形して、(ap-10)(bp-10)=100
p=1のときの解の例は、(a,b)=(20,20)
p=2のときの解の例は、(a,b)=(10,10)
p=3のときの解の例は、(a,b)=(5,10)
p=4のときの解の例は、(a,b)=(5,5)
p=5のときの解の例は、(a,b)=(4,4)
p=6のときの解の例は、(a,b)=(2,10)
p=7のときの解の例は、(a,b)=(2,5)
p=8のとき、(4a-5)(4b-5)=25 となり、(4a-5)≡(4b-5)≡3 (mod 4) だが、25は3 (mod 4)の積
に分解できないので、答えは、p=8
(2) 式を変形して、(ap-100)(bp-100)=10000
p=1〜7は、(1)の(a,b)を10倍すればよい。
p=8のときの解の例は、(a,b)=(25,25)
p=9のときの解の例は、(a,b)=(20,25)
p=10のときの解の例は、(a,b)=(20,20)
p=11のときの解の例は、(a,b)=(10,100)
p=12のときの解の例は、(a,b)=(10,50)
p=13のときの解の例は、(a,b)=(8,200)
p=14のときの解の例は、(a,b)=(10,25)
p=15のときの解の例は、(a,b)=(12,15)
p=16のとき、(4a-25)(4b-25)=625 となり、(4a-25)≡(4b-25)≡3 (mod 4) だが、625は3 (mod 4)
の積に分解できないので、答えは、p=16
(3) 式を変形して、(ap-1000)(bp-1000)=1000000
p=1〜7は、(1)の(a,b)を100倍、p=8〜15は、(2)の(a,b)を10倍すればよい。
p=16のときの解の例は、(a,b)=(125,125)
p=17のときの解の例は、(a,b)=(60,3000)
p=18のときの解の例は、(a,b)=(100,125)
p=19のときの解の例は、(a,b)=(56,875)
p=20のときの解の例は、(a,b)=(100,100)
p=21のときの解の例は、(a,b)=(50,1000)
p=22のときの解の例は、(a,b)=(50,500)
p=23のとき、(23a-1000)(23b-1000)=1000000 となり、(23a-1000)≡(23b-1000)≡12 (mod
23)
だが、1000000は12 (mod 23)の積に分解できないので、答えは、p=23
# 手計算なので間違いがあるかも知れません
GAI さんからのコメントです。(令和7年5月31日付け)
全て正解です。手計算でいけるんですか!
では、1/a+1/b=p/10^9 を満たす異なる解は何通りあるかは、brute force ではとても時間
がかかり無理と思われますが・・・。
らすかるさんからのコメントです。(令和7年6月1日付け)
もし右辺の分母が10^9のとき、23058通りで正しければ、
10^1: 20通り
10^2: 102通り
10^3: 356通り
10^4: 958通り
10^5: 2192通り
10^6: 4456通り
10^7: 8260通り
10^8: 14088通り
10^9: 23058通り
10^10: 35896通り
10^11: 53932通り
10^12: 79174通り
10^13: 112824通り
10^14: 156434通り
10^15: 215984通り
10^16: 290394通り
10^17: 384320通り
10^18: 502942通り
10^19: 646852通り
10^20: 820292通り
のようになるかと思います。
上記を求めるのに、計算方法の工夫を重ね、最終的には、(Pari/GP形式で)
f(n)=sum(k=0,n+n,sum(m=0,floor(log(10^n/2^k)/log(5)),numdiv(gcd(2^k*5^m+10^n,2^(n+n-k)*5^(n+n-m)+10^n))))
という式で求められることがわかりました。
GAI さんからのコメントです。(令和7年6月2日付け)
らすかるさんの素因数分解式を参考に、自分流でプログラムを組んでみました。
f(n)={Div=divisors(10^(2*n));X=select(i->i<=10^n,Div);Y=apply(i->10^(2*n)/i,X);
A=apply(i->i+10^n,X);B=apply(i->i+10^n,Y);}
M=[];for(n=1,#A,M=concat(M,[gcd(A[n],B[n])]));vecsum(apply(i->#divisors(i),M))
n;
20;820292(通り)
21;1038320
22;1292462
23;1590916
24;1946888
25;2359396
26;2830798
27;3393902
28;4039842
29;4775820
30;5636084
・・・・・・・・・・・・
と一瞬で求まっていくのですね。(快適!)
以下、工事中!