2021/11/30(火)NVR2021

11月27,28日に、NVR2021という研究集会を開催しました。

今回で5回めで、コロナ禍で昨年に引き続き残念ながらオンラインです。今年は、大石CRESTが最終年度で成果報告会を開催してほしいという要請があり、時期も近くて参加メンバーも多くが被るので、共同開催することにしました。CREST報告会が土曜日、NVRが日曜日でした。いろいろ重なって大変でしたが、無事終わってホッとしています。

使用したスライドを張っておきます。1つ目はdouble-doubleに関する話題、2つ目はKrawczyk法の候補者集合に関する話題です。2つ目は軽い話題で発表時間も短め(20分)です。

2021/11/12(金)kv-0.4.53

kvライブラリを0.4.53にアップデートしました。

今回は2つのbug-fixです。

1つ目は、非常に特殊な条件ですが、xがdouble-double数で、その上端が(DBL_MAX+正数)というoverflow寸前の巨大数だった場合に、sqrt(x)の計算結果の幅が広くなってしまう(double-doubleとしてふさわしい精度が出ずにdouble程度の区間幅になってしまう)というバグの修正です。実際、
#include <kv/interval.hpp>
#include <kv/dd.hpp>
#include <kv/rdd.hpp>

int main(){
        kv::interval<kv::dd> x;
        std::cout.precision(34);

        x = std::numeric_limits<kv::dd>::max();

        std::cout << x << "\n";
        std::cout << sqrt(x) << "\n";
}
を実行すると、
[1.340780792994259561100831764080293e+154,1.340780792994259672743259815885529e+154]
と計算されてしまっていました。修正後は、
[1.340780792994259672743259815885478e+154,1.340780792994259672743259815885529e+154]
ときちんと計算されるようになりました。中田真秀先生の、
のtweetで気づくことができました(kvではなくQDライブラリへの言及ではありますが)。中田先生ありがとうございました。

2つめは、自動step幅調節機能付き定積分を行う、defint_autostepで、定積分
\int_a^b f(x) dx
を計算するとき、終端値bが区間(例えば\pi)だと無限ループ(あるいは非常に長時間のループ)に陥ってしまうことがあるというバグです。この関数は、(尺取り虫のように自動的に計算精度を見ながら)積分区間を
a=t_0, t_1, t_2, \ldots, t_{n-1}, t_n=b
のように分割して積分を行っていますが、そのプロセスは2passになっていて、t_iでTaylor展開の係数を計算して適切なstep sizeを推定し、[t_i,t_{i+1}]での精度保証付きの積分値を計算し、その区間幅を見て広すぎるか狭すぎるかしたら再度step sizeを修正してもう一度精度保証付きの積分値を計算する、ということをしています。ここで、[t_i,t_{i+1}]での積分値を計算するとき、t_it_{i+1}のどちらかが(幅が広めの)区間だった場合は積分値の計算結果の区間幅が極端に大きくなり、精度を出そうとしてt_{i+1}を引き戻してstep sizeを小さくしようとします(そんなことをしても精度は上がらないのですが)。従って、bが区間だった場合には[t_i,t_{i+1}=b]と最後の一歩を計算しようとすると、精度不足でt_{i+1}が引き戻されてしまい、いつまで経ってもbに到達しない、という現象が発生していました。最後の一歩のときはstep sizeの見直しをしないように修正しました。ずっと存在していたバグですが、0.4.51以前は最後の一歩の区間幅の再計算が広めになるバグがあってたまたま発生しづらくなっていて、0.4.52で直ったことにより顕在化したようです。

この現象は、浅井 大晴さん、多田 秀介さん、宮内 洋明さんの3名の、Bessel関数J_0(x)の計算が終わらないという指摘で気付きました。Bessel関数は、\displaystyle J_n(x) = \frac{1}{\pi} \int_0^{\pi} \cos(x \sin t - nt) dt で計算しており、終端値が区間である積分をしていたので、この現象が発生しました。みなさま、ありがとうございました。
OK キャンセル 確認 その他