Computer History Museum http://www.computerhistory.org/ のサイトに、ILLIAC IV を紹介した当時(1974)のパンフレットがあります。 http://www.computerhistory.org/brochures/full_record.php?iid=doc-4372957044511 ILLIAC IV は、並列コンピュータの祖…

Implementing the Standards - including Fortran 2003 Malcolm Cohen, Numerical Algorithms Group Ltd http://www.fortran.bcs.org/2007/jubileeprog.php http://www.fortran.bcs.org/2007/jubilee/f50.pdf こんなダメダメな forallちゃんですが、取り柄も…

大体において、組み込み関数の乱数は腐っているもので、真面目な計算に使うと怖い人達に叱られてしまいます。Intel Visual Fortran コンパイラー XE 12.1.4 で軽く検証してみます。 Fortran90 で導入されたサブルーチン RANDOM_NUMBER(x) は、elemental なサ…

Intel Fortran の intrinsic function の SUM は、大きな配列の和をとる時、配列の並びで大きく値が変わる。（浮動小数点では交換則が成り立っていないが、その誤差が大きく見えてくる。）これは、倍精度変数をアキュームレータにするか、いわゆる Kahan 和…

線型方程式 Ax=b をガウス消去法で解くプログラム。Fortran90 だと式をそのままプログラムに直せばいいので楽。でもPivotting を入れるとやや見苦しくなる。並列計算機では LU分解で Crout法より Gauss 消去法の方が有利になることもあるとか。ベクトル機全…

先日、本屋で Ruby 本を立ち読みしていたのですが、実数どうしの合同の比較で、画面の出力と実行結果を一致させるために、数値を文字列に変換して文字列どうしで比較するという項目があって、面白く奇妙な発想だと思い、メモっておくことにしました。FORTRAN…

アルキメデスは、円の内接多角形と外接多角形を書いて円周の真値を挟み込む形で円周率を求める方法を後世に残しています。ここでは内接多角形から円周率を求めることにします。ところで、小学生の知識で周長のわかる場合としては、円は内接六角形と外接正方…

二進法表示で小数点以下 n+1 桁目から始まる LOG(2) のビット列表示です。計算する式 ここで {} は、小数部を意味します。 ソース・プログラム 整数に変換することでビット列の表示をするように直しました。 PROGRAM test ! LOG(2) IMPLICIT NONE INTEGER ::…

九月二日の続き。 Pi Directory by David H. Bailey下のDavid H. Bailey, "The BBP Algorithm for Pi," manuscript, Sep 2006. に基づいて、LOG(2) の任意位置での二進表示の計算というものをやってみました。 アルゴリズムを理解できたかの試し計算なので、…

Fortran90 用の Reference card が、Homepage of Michael Goerz にあります。

コンピュータによる多倍長計算を手段とする実験数学を提唱して、最近気炎を上げている D. H. Bailey という米国人がおります。そのべ氏の製作による多倍長計算パッケージを MPFUN90 といいます。High-Precision Software Directoryこのパッケージは、Fortran…

重なった領域をコピーする場合、部分配列の代入は単純なDO..LOOPとは動作が異なります。 ■実行結果 ■ソース・プログラム 配列生成子にF2003の[,]を使いました。F90の場合はそれぞれ(/,/)となります。 PROGRAM test IMPLICIT NONE INTEGER :: i, ia(10) ia = …

一次元配列a(:)に対してMINLOC(a)の類は要素1個の配列を返すが、MINLOC(a, 1)とすればスカラーで答えを受け取れる。 PROGRAM test IMPLICIT NONE REAL :: a(3) INTEGER :: k0(0), k1(1) CALL RANDOM_NUMBER(a) k0 = SHAPE( MINLOC(a, 1) ) k1 = SHAPE( MINLO…