码迷,mamicode.com
首页 > 编程语言 > 详细

欧拉计划(python) problem 64

时间:2015-02-01 16:14:25      阅读:303      评论:0      收藏:0      [点我收藏+]

标签:

Odd period square roots

Problem 64

All square roots are periodic when written as continued fractions and can be written in the form:

N = a0 +
1
  a1 +
1
    a2 +
1
      a3 + ...

For example, let us consider √23:

√23 = 4 + √23 — 4 = 4 + 
1
 = 4 + 
1
 
1
√23—4
  1 + 
√23 – 3
7

If we continue we would get the following expansion:

√23 = 4 +
1
  1 +
1
    3 +
1
      1 +
1
        8 + ...

The process can be summarised as follows:

a0 = 4,  
1
√23—4
 = 
√23+4
7
 = 1 + 
√23—3
7
a1 = 1,  
7
√23—3
 = 
7(√23+3)
14
 = 3 + 
√23—3
2
a2 = 3,  
2
√23—3
 = 
2(√23+3)
14
 = 1 + 
√23—4
7
a3 = 1,  
7
√23—4
 = 
7(√23+4)
7
 = 8 +  √23—4
a4 = 8,  
1
√23—4
 = 
√23+4
7
 = 1 + 
√23—3
7
a5 = 1,  
7
√23—3
 = 
7(√23+3)
14
 = 3 + 
√23—3
2
a6 = 3,  
2
√23—3
 = 
2(√23+3)
14
 = 1 + 
√23—4
7
a7 = 1,  
7
√23—4
 = 
7(√23+4)
7
 = 8 +  √23—4

It can be seen that the sequence is repeating. For conciseness, we use the notation √23 = [4;(1,3,1,8)], to indicate that the block (1,3,1,8) repeats indefinitely.

The first ten continued fraction representations of (irrational) square roots are:

√2=[1;(2)], period=1
√3=[1;(1,2)], period=2
√5=[2;(4)], period=1
√6=[2;(2,4)], period=2
√7=[2;(1,1,1,4)], period=4
√8=[2;(1,4)], period=2
√10=[3;(6)], period=1
√11=[3;(3,6)], period=2
√12= [3;(2,6)], period=2
√13=[3;(1,1,1,1,6)], period=5

Exactly four continued fractions, for N ≤ 13, have an odd period.

How many continued fractions for N ≤ 10000 have an odd period?


Answer:
1322
Completed on Sun, 1 Feb 2015, 06:17

Go to the thread for problem 64 in the forum.


Python code :

import math
sqrt=math.sqrt


def func(x):
    index=0
    b=[0.0,0.0]
    c=[0.0]
    intcoeff={}
    fractioncoeff={}
    a=int(sqrt(x))
    if sqrt(x)-a==0:
        return 0
    b[0]=1.0
    b[1]=-a*1.0
    c=1.0
    index+=1
    intcoeff[index]=a
    tempstr=str(b[0])+‘_‘+str(b[1])+‘_‘+str(c)
    fractioncoeff[tempstr]=index
    while 1:
        result=rationalization(b,c,x)
        a=result[‘a‘]
        b=result[‘b‘]
        c=result[‘c‘]
        tempstr=str(b[0])+‘_‘+str(b[1])+‘_‘+str(c)
        k=fractioncoeff.get(tempstr)
        if k==None:
            index+=1
            intcoeff[index]=a
            fractioncoeff[tempstr]=index
        else:
            return index-k+1
# 分母包含无理数,分子是有理数
def rationalization(b,c,x):
    result={}
    tc=b[0]*b[0]*x-b[1]*b[1]
    b[0]=int(b[0]*c)
    b[1]=int(-b[1]*c)
    a=int((b[0]*sqrt(x)+b[1])/tc)
    b[1]-=int(a*tc)
    k=gcd(gcd(b[0],abs(b[1])),tc)
    b[0]/=k
    b[1]/=k
    tc/=k
    result[‘a‘]=a
    result[‘b‘]=b
    result[‘c‘]=tc
    return result
#求最大公约数
def gcd(m,n):
    k=min(m,n)
    for i in range(k,0,-1):
        if m%i==0 and n%i==0:
            return i
N=10000
count=0
for i in range(2,N+1):
    k=int(func(i))
    if k%2==1:
        count+=1
print(count)


time: 3s

欧拉计划(python) problem 64

标签:

原文地址:http://blog.csdn.net/zhangzhengyi03539/article/details/43371483

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!