cco*_*ook 56 fortran pi fortran77
定义PI的动机是什么？
PI=4.D0*DATAN(1.D0)
在Fortran 77代码中？我理解它是如何工作的,但是,理由是什么？
Joh*_*zen 64
此样式确保在为PI分配值时使用ANY体系结构上可用的最大精度.
jas*_*son 14
因为Fortran没有内置常量PI.但是,不是手动输入数字并且可能犯错误或者没有在给定的实现上获得最大可能的精度,让库为您计算结果可以保证这些缺点都不会发生.
这些是等价的,你有时也会看到它们:
PI=DACOS(-1.D0)
PI=2.D0*DASIN(1.D0)
Jus*_*tin 13
我相信这是因为这是关于pi的最短系列.这也意味着它是最准确的.
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,而不是4*atan(1)基于Gregory-Leibniz系列.换句话说,REAL证明是:
sin ^ 2 x + cos ^ 2 x = 1 [定理]如果x = pi/4弧度,则sin ^ 2 x = cos ^ 2 x,或sin ^ 2 x = cos ^ 2 x = 1/2.
然后,sin x = cos x = 1 /(根2).tan x(sin x/cos x)= 1,atan x(1/tan x)= 1.
因此,如果atan(x)= 1,则x = pi/4,并且atan(1)= pi/4.最后,4*atan(1)= pi.)
请不要加载我的评论 - 我还是一个青少年.
这是因为这是计算pi任意精度的一种精确方法.您可以简单地继续执行该函数以获得越来越高的精度,并在任何点停止以获得近似值.
相比之下,指定pi为常量可为您提供与最初给定的精度完全相同的精度,这可能不适合高度科学或数学应用(因为Fortran经常使用).
这个问题不仅仅令人眼前一亮。为什么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)
浮点参数2：在二进制表示形式中，数字表示为0.b n b n-1 ... b 0 x 2 m。如果逆三角函数为其参数提供了最佳的数值二进制近似值，则我们不想因乘法而损失精度。为此，我们只能乘以2的幂。
? = 4 arctan(1) = 2 arcsin(1)
= 2 arccos(0)
= -4 arctan(-1) = arccos(-1) = -2 arcsin(-1)
注意：这在IEEE-754 binary64表示形式（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 binary32（最常见的形式REAL或kind=REAL32）和IEEE-754 binary128（最常见的形式kind=REAL128）
模糊实现的论点：从这一点来看，一切都取决于三角函数的实现。有时arccos和arcsin分别来自atan2和atan2作为
ACOS(x) = ATAN2(SQRT(1-x*x),1)
ASIN(x) = ATAN2(1,SQRT(1-x*x))
或更具体地从数字角度来看：
ACOS(x) = ATAN2(SQRT((1+x)*(1-x)),1)
ASIN(x) = ATAN2(1,SQRT((1+x)*(1-x)))
此外，x86指令集的atan2一部分，其他则不是。为此，我将争论以下用法：FPATAN
? = 4 arctan(1)
超过所有其他。
注意：这是一个模糊的论点。我敢肯定有人对此有更好的看法。
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任何类型的a都应表示IEEE-754浮点数。的表示REAL取决于处理器。这意味着我可以查询selected_real_kind(33, 4931)并期望获得binary128浮点数，但是我可能会kind返回一个以更高的精度表示浮点的返回值。也许是100位数，谁知道呢。在这种情况下，我上面的数字串太短了！不能肯定使用它吗？甚至那个文件可能太短了！
有趣的事实 ： sin(pi) is never zero
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