MENU
407,224

スレッドNo.2697

解が存在できない条件

ディオファントス方程式
(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を各場合について求めて下さい。

引用して返信編集・削除(未編集)

(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

# 手計算なので間違いがあるかも知れません

引用して返信編集・削除(未編集)

全て正解です。
手計算でいけるんですか!
では
1/a+1/b=p/10^9
を満たす異なる解は何通りあるかはbrute force ではとても時間がかかり無理と思われますが・・・

引用して返信編集・削除(未編集)

もし右辺の分母が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通り
のようになるかと思います。

引用して返信編集・削除(編集済: 2025年06月01日 17:40)

上記を求めるのに計算方法の工夫を重ね、最終的には
(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))))
という式で求められることがわかりました。

引用して返信編集・削除(未編集)

らすかるさんの素因数分解式を参考に、自分流でプログラムを組んでみました。
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
・・・・・・・・・・・・
と一瞬で求まっていくのですね。(快適!)

引用して返信編集・削除(編集済: 2025年06月02日 15:05)

このスレッドに返信

ロケットBBS

Page Top