基本三角関数の浮動小数点例外

Question

地理座標から地心座標に変換するサブルーチンの浮動小数点例外に問題があります。変数geo(x)は、緯度と経度のペアとしてサブルーチンに入力されます。変数xyz(x)は、構成要素（x - グリニッジ子午線および赤道; y - 経度90度および赤道90度; z - 北極）のトリプレットとして出力されます。

subroutine geo2xyz(geo,xyz) 
IMPLICIT REAL*8 (A-H,O-Z) 
dimension geo(2),xyz(3) 
rr=6367443.5 
xyz(1)=rr*sind(90.-geo(1))*cosd(geo(2)) 
xyz(2)=rr*sind(90.-geo(1))*sind(geo(2)) 
xyz(3)=rr*cosd(90.-geo(1)) 
return 
end 

私はsindとcosdは、非標準機能であることを理解、しかし、彼らは適切な変換を持つオブジェクトファイルにリンクされています。私は前にこれをテストして、それが代わりにラジアンの度を使用する他のコードのために働く：

real function sind(x) 
IMPLICIT REAL*8 (A-H,O-Z) 
sind=sin(x*3.141592653589793d0/180.0d0) 
return 
end 

real function cosd(x) 
IMPLICIT REAL*8 (A-H,O-Z) 
cosd=cos(x*3.141592653589793d0/180.0d0) 
return 
end 

私はGDBを使用してプログラムを通過しようとした、と問題を把握できませんでした。方程式のすべての変数はOKですが、xyz(1) = 1であると思われます。計算機で計算すると、そうではありません。

Program received signal SIGFPE, Arithmetic exception. 
0x0000000100001cb8 in geo2xyz (geo=..., xyz=...) at disslip.f:758 
758  xyz(1)=rr*sind(90.-geo(1))*cosd(geo(2)) 
(gdb) print xyz(1) 
$1 = 1 
(gdb) print rr 
$2 = 6367443.5 
(gdb) print geo(1) 
$3 = 50.350000000000001 
(gdb) print geo(2) 
$4 = -127.55800000000001 
(gdb) 

浮動小数点例外の原因は何ですか？私はそれがかなりシンプルだと確信しています、私はこれに非常に新しいです。

Answer 1

sind関数とcosd関数は、REAL（単精度）浮動小数点を返します。

しかし、subroutine geo2xyzにはそれについての手掛かりはありません。
IMPLICIT REAL*8では、それらの関数が倍精度浮動小数点を返すと仮定します。

したがって、使用する前にsindとcosdを実際に正しく宣言し、戻り値の型をreal * 8に変更する必要があります。

面白いことに、返された単精度値がST0レジスタで倍に昇格されたため、このプログラムは古いx86アーキテクチャで動作していた可能性があります。 SSE2とレジスタXMM0を使用している64ビットのインテルでは動作しません。そのようなアーキテクチャでは倍増することはありません。