1000桁のルート2(その3)

数列を格納する配列を整数型にしてみたところ、劇的に速くなりました。
なぜ今まで文字列型にしていたのか自分でもわかりません・・
(結果は改行加工してます)

C:\Temp>python hsqrt-1.3.py 2 1000


---------------結果---------------
1.
4142135623 7309504880 1688724209 6980785696 7187537694
8073176679 7379907324 7846210703 8850387534 3276415727
3501384623 0912297024 9248360558 5073721264 4121497099
9358314132 2266592750 5592755799 9505011527 8206057147
0109559971 6059702745 3459686201 4728517418 6408891986
0955232923 0484308714 3214508397 6260362799 5251407989
6872533965 4633180882 9640620615 2583523950 5474575028
7759961729 8355752203 3753185701 1354374603 4084988471
6038689997 0699004815 0305440277 9031645424 7823068492
9369186215 8057846311 1596668713 0130156185 6898723723
5288509264 8612494977 1542183342 0428568606 0146824720
7714358548 7415565706 9677653720 2264854470 1585880162
0758474922 6572260020 8558446652 1458398893 9443709265
9180031138 8246468157 0826301005 9485870400 3186480342
1948972782 9064104507 2636881313 7398552561 1732204024
5091227700 2269411275 7362728049 5738108967 5040183698
6836845072 5799364729 0607629969 4138047565 4823728997
1803268024 7442062926 9124859052 1810044598 4215059112
0249441341 7285314781 0580360337 1077309182 8693147101
7111168391 6581726889 4197587165 8215212822 951848847

----------------------------------

計算開始時刻:2018-06-22 00:37:51.250304
計算終了時刻:2018-06-22 02:30:03.352490
計算時間: 1:52:12.102186
         
計算時間: 7:09:50.459134

(Intel Core-i7920 @2.67GHz)

このVer 1.3をベースに開平法を実装したプログラムも作りましたので、改めて比較してみたいと思います。


# -*- coding: utf-8 -*-
#
# 平方根を求めるプログラム
# Writtend by T'z Factory Co.,Ltd
# Ver 1.0  2018/05/25
# Ver 1.1  2018/06/04  Python2系に対応するため、copy属性をlist関数に書き換え。表示を一部sys.stdout関数に書き換え
# Ver 1.2  2018/06/15  桁上がり候補を求める部分をマルチプロセス化
# Ver 1.3  2018/06/20  配列を文字列型から整数型へ変更

#==========================================
# インポートモジュール
#==========================================
try:
  import datetime, logging, sys, time
  from multiprocessing import Pool

except ImportError:
  print('必要なモジュールをインポートできませんでした。環境を確認してください。')
  sys.exit()

#==========================================
#変数定義
#==========================================

SEP =  '---------------'

#==========================================
# 上位桁から有効数字までの'0'を削除する
#==========================================
def zerosupress(a):
  i = len(a)
  j = 0
  while (i > 0 and j < 1):
    if a[i-1] == 0:
      a.pop()
    else:
      j += 1
    i -= 1
  i = len(a)
  if a[i-1] == '.':
    a.append(0)
  return(a)

#==========================================
# 配列同士の足し算
#==========================================
def add(a1, a2):
  a = list(a1)
  b = list(a2)

# 小数点以下の桁数を合わせる
  p = a.index('.') - b.index('.')
  if p > 0:
    for i in range(p):
      b.insert(0, 0)
  if p < 0:
    for i in range(-p):
      a.insert(0, 0)

# 整数部分の桁数を合わせる
  p = (len(a) - a.index('.')) - (len(b) - b.index('.'))
  if p > 0:
    for i in range(p):
      b.append(0)
  if p < 0:
    for i in range(-p):
      a.append(0)
  l_a = len(a)

# 和格納用配列初期化
  c = []
  for i in range(l_a + 1):
    c.append(0)

# キャリー初期化
  cy = 0

# 計算
  for i in range(l_a):
    if a[i] == '.':
      c[i] = '.'
    else:
      t = a[i] + b[i] + cy
      if t > 9:
        cy = 1
        t -= 10
      else:
        cy = 0
      c[i] = t

# 最後のキャリーを足し込む
  if cy == 1:
    c[i + 1] = cy
  return(c)

#==========================================
# 配列同士の掛け算
#==========================================
def multi(a1, a2):
  a = list(a1)
  b = list(a2)

# ピリオドを一旦最下位に移動(整数化)
  p_a = a.index('.')
  p_b = b.index('.')
  a.pop(p_a)
  b.pop(p_b)
  a.insert(0, '.')
  b.insert(0, '.')
  l_a = len(a)
  l_b = len(b)
  l_d = l_a - l_b
  l_p = l_a + l_b

# 桁数を合わせる
  if l_d > 0:
    for i in range(l_d):
      b.append(0)
  if l_d < 0:
    for i in range(-l_d):
      a.append(0)
  l_a = len(a)

# 積格納用配列初期化
  d = []
  for i in range(l_p):
    d.append(0)
  d.insert(0, '.')

# 計算
  for i in range(1, l_a):

# 和格納用配列初期化
    c = []
    for j in range(1, l_a):
      c.append(0)
    c.insert(0, '.')
    for k in range(b[i]):
      c = add(c, a)
    c.pop(0)
    for k in range(1, i):
      c.insert(0,0)
    c.insert(0, '.')
    d = add(d, c)
  d.pop(0)
  d.insert(p_a + p_b, '.')
  return(zerosupress(d))

#==========================================
# 配列同士の比較
#==========================================
# 戻り値
# a1 > a2 -- 2
# a1 < a2 -- 1
# a1 = a2 -- 0
# エラー --- -1
#==========================================
def compare(a1, a2):
  a = list(a1)
  b = list(a2)

# 小数点以下の桁数を合わせる
  p = a.index('.') - b.index('.')
  if p > 0:
    for i in range(p):
      b.insert(0, 0)
  if p < 0:
    for i in range(-p):
      a.insert(0, 0)

# 整数部分の桁数を合わせる
  p = (len(a) - a.index('.')) - (len(b) - b.index('.'))
  if p > 0:
    for i in range(p):
      b.append(0)
  if p < 0:
    for i in range(-p):
      a.append(0)
  l_a = len(a)
  r = -1

# 比較
  for i in range(l_a):
    j = l_a - i -1
    if a[j] > b[j]:
      r = 2
      break
    if a[j] < b[j]:
      r = 1
      break
    if a[j] == b[j]:
      r = 0
  return(r)

#==========================================
# 計算結果を表示
#==========================================
def display(a):
  if len(a) == 0:
    return(False)
  l_a = len(a)
  for i in range(l_a):
    j = l_a - i - 1
    if j > 0:
      b = a[j]
    if j == 0 and a[j] == '.':
      b = ' '
    else:
      b = a[j]
    sys.stdout.write(str(b))
  print('\r')

#==========================================
# 候補値を計算し比較結果を返す
#==========================================
# z -- 平方根
# x -- 平方根を求める値
# p -- 加算値
# k -- pをk倍する
#==========================================
def zcandidate(x, z, p, k):
  if k < 10:
    ks = ['.', k]
  else:
    ks = ['.', 0, 1]
  zc = zerosupress(add(z, multi(p, ks)))
  rc = compare(x, multi(zc, zc))
  return(rc, zc)

#==========================================
# ラッパー
#==========================================
def wrap_zc(args):
  return(zcandidate(*args))

#==========================================
# 平方根の計算
#==========================================
# x -- 平方根を求める値
# c -- 求める平方根の桁数
#==========================================
def sqrt(x, c):

# max_proc = multiprocessing.cpu_count()

  max_proc = 10
  pool = Pool(max_proc)

# 平方根の初期値
  z = ['.', 0]
  
# 加算値
  p = ['.', 1]
  m10 = ['.', 0, 1]
  i = 0
  while i < c:
    zcargs = [(x, z, p, 1), (x, z, p, 2), (x, z, p, 3), (x, z, p, 4), (x, z, p, 5), (x, z, p, 6), (x, z, p, 7), (x, z, p, 8), (x, z, p, 9), (x, z, p, 10)]
    zcr = pool.map(wrap_zc, zcargs)
    j = 0
    r = 0
    while j < 10:
# 平方根が得られた
      if zcr[j][0] == 0:
        r = 9
        return(zcr[j][1])
# 二乗値が目標値を超えた
      if zcr[j][0] == 1:
        if j == 0:
          z.insert(0, 0)
        else:
          z = zcr[j-1][1]
        r = 1
        if p[len(p)-1] == 1:
          p[len(p)-1] = 0
          p.insert(0, 1)
        else:
          p[0] = 0
          p.insert(0, 1)
        break
      j += 1

# 目標値を超えなかった
    if r == 0:
      z = zerosupress(add(z, multi(p, m10)))
    i = len(z) - 1
    if i % 10 == 0:
      print('%d桁目の計算結果です。' % i)
      print('\r')
      display(z)
      print('------------------------------')
  return(z)

#==========================================
# メイン
#==========================================
if __name__ == '__main__':
  try:
    xi = list(str(sys.argv[1]))
    c = int(sys.argv[2])
    x_len = len(xi)
    if c == 0:
      raise exception
  except:
    print('計算回数を入力してください')
    print('Usage: sqrt.py xxx count')
    sys.exit()
    
  x = []
  y = []
  p = 0
  for i in range(x_len):
    if xi[i] == '.':
      p += 1
    x.insert(0, int(xi[i]))

# ピリオドが2個以上ならエラー
  if p > 1:
    print('入力した数字が誤っています')
    sys.exit()

# ピリオドが含まれていない場合は最終桁にピリオドを追加する
  if p == 0:
    x.insert(0, '.')

# 計算開始時刻を記録
  print(SEP + SEP)
  st = datetime.datetime.now()

# 平方根を計算
  y = sqrt(x, c)

# 計算終了時刻を記録
  et = datetime.datetime.now()

# 表示
  print('\r')
  print(SEP + '結果' + SEP)
  display(y)
  print(SEP + '----' + SEP)

# 終了処理

  print('\r')
  print('\r')
  sys.stdout.write('計算開始時刻:')
  print(str(st))
  sys.stdout.write('計算終了時刻:')
  print(str(et))
  sys.stdout.write('計算時間:     ')
  print(str(et-st))

  sys.exit()

 


2018年06月24日