来訪者数が２５万人を越えている．最近は１日の来訪者が２千を超す場合もある（滅多にないが）．どこから湧いてくるのかよく分からないが，数が増えるのは単純にうれしい．
べき乗和の一般公式を求めているところだが，ゴールがやや近付いたという雰囲気が出てきた．黒木の公式とファウルハーバーの公式を架橋するというのが目標だが，黒木の公式を一般化できる可能性が出てきた．まだ，どこかに計算ミスが残っているようだが，もし，これで基本的なところが通っているとすれば，解決は時間の問題だ．そのためには少し環境整備しなくてはならない．複数の公式・恒等式を使わなくてはならないが，その中にはまだ数値的に動作確認されていないものも含まれる．その辺りは黒木氏にお願いしているところなのだが，あまり積極的ではないように見かけるので，そろそろ自分でやる必要が出てきた．まず，julia のサビ落としから始めなくてはならない．
動作確認する必要のある関数（アルゴリズム）には以下がある．
- 黒木の公式 ✓
- ファウルハーバーの公式
- べき乗和の公式 ✓
- ベルヌーイ数列の一般項
- ベルヌーイ数マトリックスの生成
- 黒木の恒等式
- べき指数の逆順変換
- べき乗和の一般公式
まず，黒木の公式から確認したいのだが，方法はあるだろうか？今の環境では変数の入った多項式を文字式として扱うことはできない．従って，できるとしても，ｋとｎを指定して値を基準値と比較するしかない．そのためには，まず第３項のべき乗和の公式を作る必要がある．これはｋごとにF_k(n)=F(n)のような多項式として表示した式だ．現在，julia 環境には，仕掛りのプロジェクトとして「階乗素因数のべき」というのが上がっている．走らせると，以下のような出力が吐き出される．
n=14^16 ゼロ=2494949494949492 n=10000000000000000 r=6.000000000000000041323199371444658202404150485
572577732799965364508923497182429
多分，これはAMSのDaily Epsilonに関係するものと思われるが，すっかり忘れてしまった．TrailingZeros という関数名があるので，巨大数末尾のゼロの個数をカウントしていたものと思われるが… その他に，primeExponent という関数も入っている．階乗数Nをｂ進数表記したときの末尾ゼロの個数をカウントしているのだろう．使っているモジュールは Printf と Primes だけだ．高精度演算が必要なので，setprecision などが使われている．データ型はBigIntとBigFloatを使っている．あとは，@printf や println だけだ．
べき乗和公式というフォルダを作って，ここから始めることにした．⇒PowSum1～5という関数セットを作った．動作している．べき乗和を単純に計算する処理を組み込んで結果が一致していることを確認した．べき乗するところはbig(n)^kのような形式で実行する必要がある．さもないと桁落ちが発生してしまう．さて，いよいよ黒木の公式だ．これは大きいのでかなり大変だろう．
いや，かなり簡単に実装できた．実際，この関数はべき乗和関数を再帰実行しているので，既成のPowerSumがそのまま使えた．これで３つの関数が出揃ったことになる．動作確認が必要なルーチンがあと，２つある．黒木の恒等式とべき指数の逆順変換だ．黒木の恒等式は片付いた．べき指数の逆順変換も通った．これで準備は整った．