当前位置： 文江博客话题详情

为什么定义 PI = 4*ATAN(1.d0)

发布于 2024-08-19 16:32:01 字数 112 浏览 5 评论 0原文

定义 PI 的动机是什么

PI=4.D0*DATAN(1.D0)

在 Fortran 77 代码中 ？我明白它是如何工作的，但是，推理是什么？

原文

What is the motivation for defining PI as

PI=4.D0*DATAN(1.D0)

within Fortran 77 code? I understand how it works, but, what is the reasoning?

如果你对这篇内容有疑问，欢迎到本站社区发帖提问 参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

扫码二维码加入Web技术交流群

发布评论

需要 登录 才能够评论， 你可以免费 注册 一个本站的账号。

评论（6）

情愿 2024-08-26 16:32:01

这种风格确保在为 PI 分配值时使用任何架构上可用的最大精度。

This style ensures that the maximum precision available on ANY architecture is used when assigning a value to PI.

回复 原文
千紇 2024-08-26 16:32:01

因为 Fortran 没有内置的 PI 常量。但是，让库为您计算结果可以保证这些缺点都不会发生，而不是手动输入数字并可能犯错误或无法在给定的实现上获得最大可能的精度。

这些是等效的，有时您也会看到它们：

PI=DACOS(-1.D0)
PI=2.D0*DASIN(1.D0)

Because Fortran does not have a built-in constant for PI. But rather than typing in the number manually and potentially making a mistake or not getting the maximum possible precision on the given implementation, letting the library calculate the result for you guarantees that neither of those downsides happen.

These are equivalent and you'll sometimes see them too:

PI=DACOS(-1.D0)
PI=2.D0*DASIN(1.D0)
回复 原文
稚气少女 2024-08-26 16:32:01

我相信这是因为这是圆周率上最短的系列。这也意味着它是最准确的。

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。

那么，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。）

请不要给我评论——我还是个青春期前的孩子。

I believe it's because this is the shortest series on pi. That also means it's the MOST ACCURATE.

The Gregory-Leibniz series (4/1 - 4/3 + 4/5 - 4/7...) equals pi.

atan(x) = x^1/1 - x^3/3 + x^5/5 - x^7/7...

So, 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...

That equals the Gregory-Leibniz series, and therefore equals pi, approximately
3.1415926535 8979323846 2643383279 5028841971 69399373510.

Another way to use atan and find pi is:

pi = 16*atan(1/5) - 4*atan(1/239), but I think that's more complicated.

I hope this helps!

(To be honest, I think the Gregory-Leibniz series was based on atan, not 4*atan(1) based on the Gregory-Leibniz series. In other words, the REAL proof is:

sin^2 x + cos^2 x = 1 [Theorem]
If x = pi/4 radians, sin^2 x = cos^2 x, or sin^2 x = cos^2 x = 1/2.

Then, sin x = cos x = 1/(root 2). tan x (sin x / cos x) = 1, atan x (1 / tan x) = 1.

So if atan(x) = 1, x = pi/4, and atan(1) = pi/4.
Finally, 4*atan(1) = pi.)

Please don't load me with comments-I'm still a pre-teen.

回复 原文
半窗疏影 2024-08-26 16:32:01

这个问题的内涵比表面上看到的要复杂得多。为什么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.bnbn-1.. .b0 x 2米。如果反三角函数为其参数提供了最佳的数值二进制近似值，我们不希望因乘法而损失精度。为此，我们应该只乘以 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 最常见的形式）

实现参数：在 intel CPU 上，atan2 是 x86 指令集的一部分，作为 FATAN而另一个反三角函数则源自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 的有趣读物： arctan 是如何实现的？, < a href="https://retrocomputing.stackexchange.com/q/14244">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) 并期望获得 binary128 浮点数，但我可能会返回一个代表精度更高的浮点数的 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

There is more to this question than meets the eye. Why 4 arctan(1)? Why not any other representation such as 3 arccos(1/2)?

This will try to find an answer by exclusion.

Mathematical intro: When using inverse trigonometric functions such as arccos, arcsin and arctan, one can easily compute π in various ways:

π = 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) = ...

There exist many other exact algebraic expressions for trigonometric values that could be used here.

Floating-point argument 1: it is well understood that a finite binary floating-point representation cannot represent all real numbers. Some examples of such numbers are 1/3, 0.97, π, sqrt(2), .... To this end, we should exclude any mathematical computation of π where the argument to the inverse trigonometric functions cannot be represented numerically. This leaves us the arguments -1,-1/2,0,1/2 and 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)

Floating-point argument 2: in the binary representation, a number is represented as 0.bnbn-1...b0 x 2m. If the inverse trigonometric function came up with the best numeric binary approximation for its argument, we do not want to lose precision by multiplication. To this end we should only multiply with powers of 2.

π = 4 arctan(1) = 2 arcsin(1)
  = 2 arccos(0)
  = -4 arctan(-1) = arccos(-1) = -2 arcsin(-1)

Note: this is visible in the IEEE-754 binary64 representation (the most common form of DOUBLE PRECISION or kind=REAL64). There we have

write(*,'(F26.20)') 4.0d0*atan(1.0d0) -> "    3.14159265358979311600"
write(*,'(F26.20)') 3.0d0*acos(0.5d0) -> "    3.14159265358979356009"

This difference is not there in IEEE-754 binary32 (the most common form of REAL or kind=REAL32) and IEEE-754 binary128 (the most common form of kind=REAL128)

Implementation argument: On intel CPU, the atan2 is part of the x86 Instruction set as FPATAN while the other inverse trigonometric function are derived from atan2. A potential derivation could be:

          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)))

This can be seen in the assembly code of these instructions (See here). To this end I would argue the usage of:

π = 4 arctan(1)

Note: this is a fuzzy argument. I'm certain there are people with better opinions on this.
Interesting reads on FPATAN: How is arctan implemented?, x87 trigonometric instructions

The Fortran argument: why should we approximate π as :

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)

and not :

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

The answer lays in the Fortran standard. The standard never states that a REAL of any kind should represent an IEEE-754 floating point number. The representation of REAL is processor dependent. This implies that I could inquire selected_real_kind(33, 4931) and expect to obtain a binary128 floating-point number, but I might get a kind returned that represents a floating-point with much higher accuracy. Maybe 100 digits, who knows. In this case, my above string of numbers is to short! One cannot use this just to be sure? Even that file could be too short!

Interesting fact : 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"

which is understood as:

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
回复 原文
坦然微笑 2024-08-26 16:32:01

这是因为这是一种以任意精度计算 pi 的精确方法。您可以简单地继续执行该函数以获得越来越高的精度，并在任意点停止以获得近似值。

相比之下，将 pi 指定为常量可为您提供与最初给出的精确度完全相同的精度，这可能不适合高度科学或数学应用程序（正如 Fortran 经常使用的那样）。

It's because this is an exact way to compute pi to arbitrary precision. You can simply continue executing the function to get greater and greater precision and stop at any point to have an approximation.

By contrast, specifying pi as a constant provides you with exactly as much precision as was originally given, which may not be appropriate for highly scientific or mathematical applications (as Fortran is frequently used with).

回复 原文
成熟稳重的好男人 2024-08-26 16:32:01

这听起来很像编译器错误的解决方法。或者这个特定的程序可能依赖于该身份的精确性，因此程序员做出了保证。

That sounds an awful lot like a work-around for a compiler bug. Or it could be that this particular program depends on that identity being exact, and so the programmer made it guaranteed.

回复 原文
~没有更多了~

关于作者

￠好甜

暂无简介

0 文章
0 评论
22 人气
发私信

相关话题

更多

热门标签

操作系统 程序设计 IT运维 Linux系统管理 JavaScript 服务器应用 solaris C/C++ PHP Shell BSD Vue.js aix Oracle Python HTML 系统管理 HTML5 CSS 前端
更多

推荐作者

1CH1MKgiKxn9p

文章 0 评论 0

ゞ记忆︶ㄣ

文章 0 评论 0

JackDx

文章 0 评论 0

信远

文章 0 评论 0

yaoduoduo1995

文章 0 评论 0

霞映澄塘

文章 0 评论 0

更多

友情链接

文江博客
    我们使用 Cookies 和其他技术来定制您的体验包括您的登录状态等。通过阅读我们的 隐私政策 了解更多相关信息。 单击 接受 或继续使用网站，即表示您同意使用 Cookies 和您的相关数据。
    原文