X



トップページ数学
1002コメント415KB

【数セミ】エレガントな解答をもとむ3【2018.10】

■ このスレッドは過去ログ倉庫に格納されています
0774132人目の素数さん
垢版 |
2020/08/15(土) 03:00:41.98ID:Bp8epmHH
>>773
検証プログラム作成時

1より大きい無理数αに対し正規連分数の打ち切り近似有理数列ri(α)を
r0(α)=[α]、
r(i+1)(α)=[α]+1/ri(1/(α-[α]))
で定める

補題
任意の奇数iに対しri(α)はαの上からの近似であり、分母がri(α)の分母以下である任意の上からの近似の中で最良である

2/(2^(1/n)-1)の天井と2n/(log2)の床が一致しないものが存在するならば、その最小を与えるものはある奇数iについてri(2/log2)の分母として得られる自然数に限られる
主張
奇数iに対し
ri(2/log2)=p/q、
r(i+1)(2/log2)=r/s、
r(i+2)(2/log2)=t/u
とおく
1) qが反例を与えるには(p/q-2/log2)q^2<log2/6が必要である。特に
(p/q-t/u)q^2<log2/6‥@
が必要である
2)qが反例を与えるには(p/q-2/log2)q^2<log2/6-(log2)^3/360/q^2なら十分である。特に
(p/q-r/s+1/q^4)q^2<log2/6‥A
なら十分である
定理
最小の反例はr39(2/log2)の分母である
n= 777451915729368
である
∵) r1(2/log2)=3/1であるがn=1は明らかに反例を与えない
i:3〜37である奇数に対し主張の@を満たすものはない
i=39のとき主張のAは満たされる
0775132人目の素数さん
垢版 |
2020/08/15(土) 03:01:14.23ID:Bp8epmHH
import Data.Ratio

-- 一般化連分数から有理数近似列を作成
fracs bs cs = id
$ map (head.fst)
$ iterate (\([[p,q],[p',q']],n) -> ([[(bs!!n)*p+(cs!!n)*p',(bs!!n)*q+(cs!!n)*q'],[p,q]],n+1))
$ ([[bs!!0,1],[1,0]],1)

rats bs cs = id
$ map (\[x,y] -> x%y)
$ fracs bs cs

-- 1/(log2) の一般化連分数の作成
bs = [1..]
cs = 0:(concat [[k^2,k^2] | k<-[1..]])

-- 有理数近似列から正規連分数の作成
makeRbs ratapxs = let
rcfsMain rs = let
apx1 = truncate $ rs !! 0
apx2 = truncate $ rs !! 1
newrs = map (recip.(+(-(fromInteger $ apx1)))) rs
in if apx1 == apx2
then (apx1:(rcfsMain newrs))
else rcfsMain $ tail rs
in rcfsMain ratapxs

-- 1/(log2) の正規連分数の作成
rbs = makeRbs $ drop 2 $ rats bs cs
rcs = repeat 1

-- 1/(log2) の正規連分数から有理数近似列を作成
rrats = map (*2) $ rats rbs rcs

-- 反例を与えるための必要条件 ind1<(log2)/6 と
-- 反例を与えるための十分条件 ind2<(log2)/6 の
-- 検証
tests = [ (k,ind1,ind2,n) |
k<-[1,3..],
let x = rrats !! k,
let y = rrats !! (k+1),
let z = rrats !! (k+2),
let n = denominator x,
let rn = fromInteger $ n,
let ind1 = fromRational $ (x-z)*rn^2,
let ind2 = fromRational $ (x-y+1/rn^4)*rn^2
]

main = do
mapM_ print $ take 20 $ tests
0776132人目の素数さん
垢版 |
2020/08/15(土) 03:01:34.76ID:Bp8epmHH
(1,0.1111111111111111,1.1428571428571428,1)
(3,0.28125,0.3074276462254604,9)
(5,0.2157303370786517,0.37955514726065,96)
(7,0.5604534005037783,0.778657718283889,445)
(9,0.1838388515860153,0.2907375803451737,794)
(11,1.0174722211315717,1.2252482805588467,4319)
(13,0.1773551664467339,0.19461297624934124,25469)
(15,0.6869730361965083,1.0463959311459805,287209)
(17,0.13674228263475482,0.18821662431664307,418079)
(19,0.7335370174413465,1.1584026182264542,3057423)
(21,1.0170018170613446,1.2245661504073762,8336111)
(23,0.19155865628150245,0.21960325986464296,24590254)
(25,1.4623220085717246,1.5494690343620314,385107953)
(27,0.2731271551989445,0.28278157095101175,6847196937)
(29,0.8922736690121698,1.6110001975243693,200557046245)
(31,1.1787446582422778,1.287114621424387,449541554817)
(33,0.210010870755119,0.23465044264678134,2669611999270)
(35,0.2248455561686892,0.40858063754648183,25423560120208)
(37,0.7271910471480607,0.7842187311800556,113071214541301)
(39,3.29262416686285e-2,3.524736840949706e-2,777451915729368)
0777132人目の素数さん
垢版 |
2020/08/16(日) 02:41:58.19ID:ocOa8mpd
>>774
コード間違ってた
1/log2正規連分数の打ち切り近似を2倍しても2/log2の正規連分数の打ち切り近似にはならん
コード直したらn=777451915729368はi=35での上から評価だった

import Data.Ratio

-- 一般化連分数から有理数近似列を作成
fracs bs cs = id
$ map (head.fst)
$ iterate (\([[p,q],[p',q']],n) -> ([[(bs!!n)*p+(cs!!n)*p',(bs!!n)*q+(cs!!n)*q'],[p,q]],n+1))
$ ([[bs!!0,1],[1,0]],1)

rats bs cs = id
$ map (\[x,y] -> x%y)
$ fracs bs cs

-- 2/(log2) の一般化連分数の作成
bs = 2:[2..]
cs = [0,2,1]++((concat [[k^2,k^2] | k<-[2..]]))

-- 有理数近似列から正規連分数の作成
makeRbs ratapxs = let
makeRbsMain rs = let
apx1 = truncate $ rs !! 0
apx2 = truncate $ rs !! 1
newrs = map (recip.(+(-(fromInteger $ apx1)))) rs
in if apx1 == apx2
then (apx1:(makeRbsMain newrs))
else makeRbsMain $ tail rs
in makeRbsMain ratapxs

-- 2/(log2) の正規連分数の作成
rbs = makeRbs $ drop 10 $ rats bs cs
rcs = repeat 1

-- 1/(log2) の正規連分数から有理数近似列を作成
rrats = rats rbs rcs
--rrats = map (*2) $ rats rbs rcs

-- 反例を与えるための必要条件 ind1<(log2)/6 と
-- 反例を与えるための十分条件 ind2<(log2)/6 の
-- 検証
tests = [ (k,ind1,ind2,n) |
k<-[1,3..],
let x = rrats !! k,
let y = rrats !! (k+1),
let z = rrats !! (k+2),
let n = denominator x,
let rn = fromInteger $ n,
let ind1 = fromRational $ (x-z)*rn^2,
let ind2 = fromRational $ (x-y+1/rn^4)*rn^2
]

main = do
mapM_ print $ take 18 $ tests
0778132人目の素数さん
垢版 |
2020/08/16(日) 02:42:30.06ID:ocOa8mpd
(1,0.1111111111111111,1.125,1)
(3,0.2571428571428571,0.3584995251661918,9)
(5,0.3645833333333333,0.5745868183338909,35)
(7,0.24181360201511334,0.2751801401822668,96)
(9,0.21822607876241706,0.22524981315418816,794)
(11,0.18275732576857484,0.19461297624934124,25469)
(13,0.15236871728780266,0.1584026182320683,418079)
(15,0.4463335352290383,0.8061415376018488,10975455)
(17,0.18006122736094865,0.21960325986464296,24590254)
(19,0.5385099766993894,0.5494690343620313,136566069)
(21,0.24775175562084492,0.3293483999092869,6847196937)
(23,0.36333931797114377,0.5706953927377529,27637329632)
(25,0.5300370170033072,0.6110001975243692,76064791959)
(27,0.6450755494730714,0.8217786553560967,574033809103)
(29,0.210010870755119,0.23465044264678134,2669611999270)
(31,0.2533250635740805,0.2900654933452463,25423560120208)
(33,0.5163505675436897,0.6960546827950139,200718868962394)
(35,3.29262416686285e-2,3.404729100026536e-2,777451915729368)
■ このスレッドは過去ログ倉庫に格納されています

ニューススポーツなんでも実況