sample\basic\sample_bigdec.hsp » Plain Format
;============================================================
; sample_bigdec.hsp — 任意精度小数 (BigDecimal) サンプル
;
; 1. 1/3 を小数点以下 30 桁まで
; 2. √2 を 50 桁まで
; 3. Machin's formula で円周率 π を 50 桁まで
; π = 16 * arctan(1/5) - 4 * arctan(1/239)
; arctan(1/x) = 1/x - 1/(3x^3) + 1/(5x^5) - ...
;============================================================
#include "iron_bigdec.hsp"
screen 0, 800, 600
title "BigDecimal sample — Machin's formula"
font "MS Gothic", 14
objsize 200, 28
; --- 1) 1/3 を 30 桁 ----------------------------------------
mes "--- 1/3 (30 decimals) ---"
a = bigdec("1")
b = bigdec("3")
q = bigdec_div(a, b, 30, BIGDEC_HALF_UP)
mes bigdec_str(q)
bigdec_release a
bigdec_release b
bigdec_release q
mes ""
; --- 2) √2 を 50 桁 -----------------------------------------
mes "--- sqrt(2) (50 decimals) ---"
two = bigdec("2")
r = bigdec_sqrt(two, 50)
mes bigdec_str(r)
bigdec_release two
bigdec_release r
mes ""
; --- 3) Machin's formula で π を 50 桁 ----------------------
mes "--- pi via Machin's formula (50 decimals) ---"
PREC = 60 ; 計算は少し余裕を持って 60 桁、最後に 50 桁へ丸める
sixteen = bigdec("16")
four = bigdec("4")
atan_1_5 = calc_atan_reciprocal(5, PREC)
atan_1_239 = calc_atan_reciprocal(239, PREC)
t1 = bigdec_mul(sixteen, atan_1_5)
t2 = bigdec_mul(four, atan_1_239)
pi = bigdec_sub(t1, t2)
pi50 = bigdec_round(pi, 50, BIGDEC_HALF_UP)
mes bigdec_str(pi50)
mes "Reference:"
mes "3.14159265358979323846264338327950288419716939937510"
bigdec_release sixteen
bigdec_release four
bigdec_release atan_1_5
bigdec_release atan_1_239
bigdec_release t1
bigdec_release t2
bigdec_release pi
bigdec_release pi50
stop
;============================================================
; arctan(1/x) = 1/x - 1/(3x^3) + 1/(5x^5) - ...
; x: 整数, prec: 途中計算桁数
; 戻り値: 新しい BigDecimal ハンドル (呼び出し側で解放)
;============================================================
#deffunc calc_atan_reciprocal int x, int prec, local one, local xb, local x2, local term_num, \
local sum_, local sign_, local k, local denom, local cur, \
local prev_sum, local dcmp, local next_sum, local mult, \
local div_result, local str_buf
one = bigdec("1")
xb = bigdec(str(x))
x2 = bigdec_mul(xb, xb) ; x^2
; 初項: 1/x (精度 prec)
term_num = bigdec_div(one, xb, prec, BIGDEC_HALF_EVEN)
; sum_ を初項で初期化
sum_ = bigdec_copy(term_num)
k = 1 ; 次の分母係数は (2k+1)
repeat
; next term: 前の分子を x^2 で割って、符号反転して (2k+1) で割る
div_result = bigdec_div(term_num, x2, prec, BIGDEC_HALF_EVEN)
bigdec_release term_num
term_num = div_result
; 符号反転 & (2k+1) で割る
denom_str = "" + (2*k + 1)
denom = bigdec(denom_str)
cur = bigdec_div(term_num, denom, prec, BIGDEC_HALF_EVEN)
bigdec_release denom
; 加算 or 減算
if k \ 2 == 1 {
next_sum = bigdec_sub(sum_, cur)
} else {
next_sum = bigdec_add(sum_, cur)
}
bigdec_release cur
; 収束判定: |term_num / (2k+1)| が 10^(-prec) 未満なら打ち切り
; ここでは term_num の大きさで判定
; prev_sum と next_sum が同じ十進文字列になったら収束
prev_sum = sum_
sum_ = next_sum
if bigdec_precision(term_num) == 1 {
; term_num が 0 になったら終了
; ただし scale を考慮する必要があるので文字列で判定
sdim tbuf, 4096
bigdec_to_plain_str term_num, tbuf, 4096
; 全部 0 か . しか無ければ終了
ok = 1
ln = strlen(tbuf)
idx = 0
repeat ln
c = peek(tbuf, idx)
if c != '0' && c != '.' && c != '-' : ok = 0 : break
idx++
loop
if ok : bigdec_release prev_sum : break
}
; |prev_sum - next_sum| が十分小さい (scale=prec+1 の 10^(-prec) 未満) なら終了
; 単純化: 差のハンドル precision と scale を比べる
diff = bigdec_sub(prev_sum, sum_)
diff_abs = bigdec_abs(diff)
bigdec_release diff
; diff_abs == 0 なら収束
tiny = bigdec("0")
cmp = bigdec_cmp(diff_abs, tiny)
bigdec_release tiny
bigdec_release diff_abs
bigdec_release prev_sum
if cmp == 0 : break
k++
if k > 5000 : break ; 安全装置
loop
bigdec_release one
bigdec_release xb
bigdec_release x2
bigdec_release term_num
return sum_