将PI定义为什么的动机是什么？
PI=4.D0*DATAN(1.D0)
在Fortran 77代码内部？我知道它是如何工作的，但是推理是什么？
6个回答
68
这种样式确保在将值分配给PI时，使用了任何架构可用的最大精度。
- John Gietzen
8
11警告：如果你采用这种方法，请注意--并非所有编译器和底层数学库在处理浮点数的三角函数时返回的结果都是相同的，尤其是在这些函数的关键点上。最近在comp.lang.fortran上进行了长时间的讨论，那里有大部分Fortran专家。他们的结论是--指定常量为pi = 3.14159...（足够精确的位数再加上一些安全位数）。 - High Performance Mark
7高性能标记(High Performance Mark)：如果您能提供您提到的comp.lang.fortran线程的链接，那将非常好！ - jvriesem
2如果有人试图计算，例如sin(31416)，正确的数学值应该是sin(31416减去(10000乘以数学pi))，而不是sin(31416减去10000次最接近数学pi的浮点值)。 - supercat
1@jvriesem：个人而言，我怀疑在 x 超出 +/-pi/2 范围时调用 sin(x) 的情况非常少，而不是满足于最接近 x 的值的正弦函数。尽管如此，Java 语言运行时会添加额外的代码到其三角函数实现中，当要计算 e.g. 2π 附近的 sin(x) 时，它将首先减去一个精确可表示的值，该值接近于 2π 并减去该值与数学 2π 的差。 - supercat
1@jvriesem：如果想象一个浮点系统的精度为五个小数位，被要求计算sin（3.1416），它可能会（注意到π约为3.1415 + 0.000092653 + 0.00000000058979），从该值中减去3.1415（得到0.0001），然后再从中减去0.000092653（得到0.000007347），再从中减去0.00000000058979（得到0.0000073464），然后对其取正弦，最终结果为-0.0000073464。 - supercat
显示剩余3条评论
16
因为Fortran没有内置常量PI。但是，与其手动输入数字并可能出现错误或不能在给定实现上获得最大可能精度，让库为您计算结果可以保证不会发生这些缺点。
以下代码也是等价的，有时你也会看到它们：
PI=DACOS(-1.D0)
PI=2.D0*DASIN(1.D0)
- jason
1
3请注意，在Fortran 77及以后的版本中，通用名称ACOS和ASIN比特定名称DACOS和DASIN更受青睐。 - Vladimir F Героям слава
14
我认为这是有关圆周率的最短级数，因此也是最准确的。
Gregory-Leibniz级数（4/1 - 4/3 + 4/5 - 4/7...）等于pi。
atan(x) = x^1/1 - x^3/3 + x^5/5 - x^7/7...
所以，atan(1) = 1/1 - 1/3 + 1/5 - 1/7 + 1/9... 4 * atan(1) = 4/1 - 4/3 + 4/5 - 4/7 + 4/9...
那就等于Gregory-Leibniz级数，因此约等于pi 3.1415926535 8979323846 2643383279 5028841971 69399373510。
使用atan另一种找到pi的方法是：
pi = 16*atan(1/5) - 4*atan(1/239)，但我认为这更复杂。
希望能对你有所帮助！
(老实说，我认为Gregory-Leibniz级数是基于atan而不是基于Gregory-Leibniz级数基于4*atan(1)。换句话说，真正的证明是：
sin^2 x + cos^2 x = 1 [定理]
如果x = pi/4 弧度，则sin^2 x = cos^2 x，或sin^2 x = cos^2 x = 1/2。
因此如果atan(x) = 1，则x = pi/4，atan(1) = pi/4。 最后，4*atan(1) = pi。)
请不要给我留下评论-我仍然是未成年人。
- Justin
2
我不明白如何从 atan x (1/tan x) = 1 推导出 atan(x) = 1。 - lindsaymacvean
2使用公式
pi = 16*atan(1/5) - 4*atan(1/239) 并不是一个数值精确的好选择。无论是 1/5 还是 1/239 都无法用浮点数精确表示。因此，由于这些浮点数的近似，atan 已经引入了误差。 - kvantour11
这个问题比表面上看到的更加复杂。为什么是4 arctan(1)？为什么不是其他表示，例如3 arccos(1/2)?
通过排除法来寻找答案。
数学介绍：当使用反三角函数如arccos, arcsin和arctan时，可以通过多种方式轻松计算π:
π = 4 arctan(1) = arccos(-1) = 2 arcsin(1) = 3 arccos(1/2) = 6 arcsin(1/2)
= 3 arcsin(sqrt(3)/2) = 4 arcsin(sqrt(2)/2) = ...
这里还有许多其他三角函数值的精确代数表达式，可以在这里使用。
浮点参数1：众所周知，有限的二进制浮点表示不能表示所有实数。一些这样的数字例子是1/3, 0.97, π, sqrt(2), ...。因此，我们应该排除任何计算π的数学计算，其中反三角函数的参数无法以数字方式表示。这使得我们只剩下参数-1,-1/2,0,1/2和1。
π = 4 arctan(1) = 2 arcsin(1)
= 3 arccos(1/2) = 6 arcsin(1/2)
= 2 arccos(0)
= 3/2 arccos(-1/2) = -6 arcsin(-1/2)
= -4 arctan(-1) = arccos(-1) = -2 arcsin(-1)
π = 4 arctan(1) = 2 arcsin(1)
= 2 arccos(0)
= -4 arctan(-1) = arccos(-1) = -2 arcsin(-1)
注意: 这在 IEEE-754二进制64位 表示中可见（这是最常见的DOUBLE PRECISION或kind=REAL64形式）。 在那里我们有
write(*,'(F26.20)') 4.0d0*atan(1.0d0) -> " 3.14159265358979311600"
write(*,'(F26.20)') 3.0d0*acos(0.5d0) -> " 3.14159265358979356009"
这种差异在IEEE-754二进制32（最常见的
REAL或kind=REAL32形式）和IEEE-754二进制128（最常见的kind=REAL128形式）中不存在。
实现参数：在英特尔CPU上，atan2是x86指令集中的一部分，作为FPATAN，而其他反三角函数则来源于atan2。可能的推导如下： mathematically numerically
ACOS(x) = ATAN2(SQRT(1-x*x),1) = ATAN2(SQRT((1+x)*(1-x)),1)
ASIN(x) = ATAN2(1,SQRT(1-x*x)) = ATAN2(1,SQRT((1+x)*(1-x)))
这可以在这些指令的汇编代码中看到（请参见此处）。为此，我认为使用以下内容是合适的：
π = 4 arctan(1)
注意: 这是一个模糊的论点。我相信有更好的意见。
FPATAN 有趣的阅读资料: 如何实现反正切函数?, x87三角函数指令
Fortran 论点: 为什么我们应该将 π 近似为:
integer, parameter :: sp = selected_real_kind(6, 37)
integer, parameter :: dp = selected_real_kind(15, 307)
integer, parameter :: qp = selected_real_kind(33, 4931)
real(kind=sp), parameter :: pi_sp = 4.0_sp*atan2(1.0_sp,1.0_sp)
real(kind=dp), parameter :: pi_dp = 4.0_dp*atan2(1.0_dp,1.0_dp)
real(kind=qp), parameter :: pi_qp = 4.0_qp*atan2(1.0_qp,1.0_qp)
并不是：
而不是：
real(kind=sp), parameter :: pi_sp = 3.14159265358979323846264338327950288_sp
real(kind=dp), parameter :: pi_dp = 3.14159265358979323846264338327950288_dp
real(kind=qp), parameter :: pi_qp = 3.14159265358979323846264338327950288_qp
答案在Fortran标准中。标准从未规定任何一种
REAL应该表示IEEE-754浮点数。 REAL的表示是处理器相关的。这意味着我可以查询selected_real_kind(33, 4931)并期望获得一个二进制128位浮点数，但我可能会得到一个代表更高精度浮点的kind。也许是100位数字，谁知道呢。在这种情况下，我的上面的数字串就太短了！难道我们不能使用这个来确保吗？即使那个文件也可能太短！有趣的事实：
sin(pi)从不为零。write(*,'(F17.11)') sin(pi_sp) => " -0.00000008742"
write(*,'(F26.20)') sin(pi_dp) => " 0.00000000000000012246"
write(*,'(F44.38)') sin(pi_qp) => " 0.00000000000000000000000000000000008672"
这被理解为：
pi = 4 ATAN2(1,1) = π + δ
SIN(pi) = SIN(pi - π) = SIN(δ) ≈ δ
program print_pi
! use iso_fortran_env, sp=>real32, dp=>real64, qp=>real128
integer, parameter :: sp = selected_real_kind(6, 37)
integer, parameter :: dp = selected_real_kind(15, 307)
integer, parameter :: qp = selected_real_kind(33, 4931)
real(kind=sp), parameter :: pi_sp = 3.14159265358979323846264338327950288_sp
real(kind=dp), parameter :: pi_dp = 3.14159265358979323846264338327950288_dp
real(kind=qp), parameter :: pi_qp = 3.14159265358979323846264338327950288_qp
write(*,'("SP "A17)') "3.14159265358..."
write(*,'(F17.11)') pi_sp
write(*,'(F17.11)') acos(-1.0_sp)
write(*,'(F17.11)') 2.0_sp*asin( 1.0_sp)
write(*,'(F17.11)') 4.0_sp*atan2(1.0_sp,1.0_sp)
write(*,'(F17.11)') 3.0_sp*acos(0.5_sp)
write(*,'(F17.11)') 6.0_sp*asin(0.5_sp)
write(*,'("DP "A26)') "3.14159265358979323846..."
write(*,'(F26.20)') pi_dp
write(*,'(F26.20)') acos(-1.0_dp)
write(*,'(F26.20)') 2.0_dp*asin( 1.0_dp)
write(*,'(F26.20)') 4.0_dp*atan2(1.0_dp,1.0_dp)
write(*,'(F26.20)') 3.0_dp*acos(0.5_dp)
write(*,'(F26.20)') 6.0_dp*asin(0.5_dp)
write(*,'("QP "A44)') "3.14159265358979323846264338327950288419..."
write(*,'(F44.38)') pi_qp
write(*,'(F44.38)') acos(-1.0_qp)
write(*,'(F44.38)') 2.0_qp*asin( 1.0_qp)
write(*,'(F44.38)') 4.0_qp*atan2(1.0_qp,1.0_qp)
write(*,'(F44.38)') 3.0_qp*acos(0.5_qp)
write(*,'(F44.38)') 6.0_qp*asin(0.5_qp)
write(*,'(F17.11)') sin(pi_sp)
write(*,'(F26.20)') sin(pi_dp)
write(*,'(F44.38)') sin(pi_qp)
end program print_pi
- kvantour
3
除了参数sin(pi) /= 0之外，我一直不愿意假设acos(-1d0)与4 * atan(1d0)一样准确，尽管这取决于实现。 - tim18
我相信有一种经过验证的算法可以为sin()生成正确舍入的大参数值，尽管它需要相当多的额外执行时间。我没有看到参考资料。一些IBM库默认为您提供此功能。任何高质量的实现都将涉及一些额外精度的模拟，这避免了在实践中可能有用的参数（例如+-20 Pi）上的重大精度降级，或者内部m387实现从返回准确值到返回参数发生突变的（较小的）点。 - tim18
顺便提一下，在m387固件中的内部Pi常数被宣传为具有66位精度。实际上，在实现中，这并不是什么特别的；正确舍入的66位精度值具有2个低阶零位。高质量的C编译器将为您提供具有21位数字的long double常数的此值。 - tim18
9
这是一种计算任意精度 pi 的精确方法。您可以继续执行该函数以获得更高的精度，并在任何时候停止以获得近似值。
相比之下，将pi指定为常量只会提供与最初给定的精度完全相同的精度，这可能不适用于高度科学或数学应用（如经常使用的Fortran）。
- John Feminella
-9
那听起来非常像是对编译器漏洞的解决方法。或者这个特定程序依赖于该身份信息的准确性，因此程序员做了保证。
- Andrew McGregor
1
3这实际上是一种非常常见的设置PI值的方式 - 不仅在Fortran中，而且在其他语言中也是如此。 （请参见上面的评论。） - jvriesem
网页内容由stack overflow 提供, 点击上面的可以查看英文原文，
原文链接
原文链接
DATAN()函数：根据参数，ATAN()函数会自动别名为其相应的单精度和双精度版本。 - jvriesem