由Python运算π的值深入Python中科学计算的实现_python

π是一个无数人追随的真正的神奇数字。我不是很清楚一个永远重复的无理数的迷人之处。在我看来,我乐于计算π,也就是计算π的值。因为π是一个无理数,它是无限的。这就意味着任何对π的计算都仅仅是个近似值。如果你计算100位,我可以计算101位并且更精确。迄今为止,有些人已经选拔出超级计算机来试图计算最精确的π。一些极值包括 计算π的5亿位。你甚至能从网上找到包含 π的一百亿位的文本文件(注意啦!下载这个文件可能得花一会儿时间,并且没法用你平时使用的记事本应用程序打开。)。对于我而言,如何用几行简单的Python来计算π才是我的兴趣所在。
你总是可以 使用 math.pi 变量的 。它被 包含在 标准库中, 在你试图自己 计算它之前,你应该去使用它 。 事实上 , 我们将 用它来计算 精度 。作为 开始, 让我们看 一个 非常直截了当的 计算Pi的 方法 。像往常一样,我将使用Python 2.7,同样的想法和代码可能应用于不同的版本。我们将要使用的大部分算法来自Pi WikiPedia page并加以实现。让我们看看下面的代码:
 

importsys
importmath

defmain(argv):

  iflen(argv) !=1:
    sys.exit('Usage: calc_pi.py <n>')

  print'\nComputing Pi v.01\n'

  a=1.0
  b=1.0/math.sqrt(2)
  t=1.0/4.0
  p=1.0

  foriinrange(int(sys.argv[1])):
    at=(a+b)/2
    bt=math.sqrt(a*b)
    tt=t-p*(a-at)**2
    pt=2*p

    a=at;b=bt;t=tt;p=pt

  my_pi=(a+b)**2/(4*t)
  accuracy=100*(math.pi-my_pi)/my_pi

  print"Pi is approximately: "+str(my_pi)
  print"Accuracy with math.pi: "+str(accuracy)

if__name__=="__main__":
  main(sys.argv[1:])

这是个非常简单的脚本,你可以下载,运行,修改,和随意分享给别人。你能够看到类似下面的输出结果: 

 你会发现,尽管 n 大于4 ,我们逼近 Pi 精度却没有多大的提升。 我们可以猜到即使 n的值更大,同样的事情(pi的逼近精度没有提升)依旧会发生。幸运的是,有不止一种方法来揭开这个谜。使用 Python Decimal (十进制)库,我们可以就可以得到更高精度的值来逼近Pi。让我们来看看库函数是如何使用的。这个简化的版本,可以得到多于11位的数字 通常情况小Python 浮点数给出的精度。下面是Python Decimal 库中的一个例子 :

wpid-python_decimal_example-2013-05-28-12-54.png

看到这些数字。不对! 我们输入的仅是 3.14,为什么我们得到了一些垃圾(junk)? 这是内存垃圾(memory junk)。 简单点说,Python给你你想要的十进制数,再加上一点点额外的值。 只要精度小于垃圾数,它不会影响任何计算。通过设置getcontext().prec 你可以的到你想要的位数 。我们试试。

看到这些数字。不对! 我们输入的仅是 3.14,为什么我们得到了一些垃圾(junk)? 这是内存垃圾(memory junk)。 简单点说,Python给你你想要的十进制数,再加上一点点额外的值。 只要精度小于垃圾数,它不会影响任何计算。通过设置getcontext().prec 你可以的到你想要的位数 。我们试试。

很好。 现在让我们 试着用这个 来 看看我们是否能 与我们以前的 代码 有更好的 逼近 。 现在, 我通常 是反对 使用“ from library import * ” , 但在这种情况下, 它会 使代码 看起来更漂亮 。
 

importsys
importmath
fromdecimalimport*

defmain(argv):

  iflen(argv) !=1:
    sys.exit('Usage: calc_pi.py <n>')

  print'\nComputing Pi v.01\n'

  a=Decimal(1.0)
  b=Decimal(1.0/math.sqrt(2))
  t=Decimal(1.0)/Decimal(4.0)
  p=Decimal(1.0)

  foriinrange(int(sys.argv[1])):
    at=Decimal((a+b)/2)
    bt=Decimal(math.sqrt(a*b))
    tt=Decimal(t-p*(a-at)**2)
    pt=Decimal(2*p)

    a=at;b=bt;t=tt;p=pt

  my_pi=(a+b)**2/(4*t)
  accuracy=100*(Decimal(math.pi)-my_pi)/my_pi

  print"Pi is approximately: "+str(my_pi)
  print"Accuracy with math.pi: "+str(accuracy)

if__name__=="__main__":
  main(sys.argv[1:])

 
输出结果: 

 好了。我们更准确了,但看起来似乎有一些舍入。从n = 100和n = 1000,我们有相同的精度。现在怎么办?好吧,现在我们来求助于公式。到目前为止,我们计算Pi的方式是通过对几部分加在一起。我从DAN 的关于Calculating Pi 的文章中发现一些代码。他建议我们用以下3个公式:

    Bailey–Borwein–Plouffe 公式
   Bellard的公式
    Chudnovsky 算法

让我们从Bailey–Borwein–Plouffe 公式开始。它看起来是这个样子: 

 在代码中我们可以这样编写它:
 

import sys
import math
from decimal import *

def bbp(n):
  pi=Decimal(0)
  k=0
  while k < n:
    pi+=(Decimal(1)/(16**k))*((Decimal(4)/(8*k+1))-(Decimal(2)/(8*k+4))-(Decimal(1)/(8*k+5))-(Decimal(1)/(8*k+6)))
    k+=1
  return pi

def main(argv):

    if len(argv) !=2:
    sys.exit('Usage: BaileyBorweinPlouffe.py <prec> <n>')

  getcontext().prec=(int(sys.argv[1]))
  my_pi=bbp(int(sys.argv[2]))
  accuracy=100*(Decimal(math.pi)-my_pi)/my_pi

  print"Pi is approximately "+str(my_pi)
  print"Accuracy with math.pi: "+str(accuracy)

if __name__=="__main__":
  main(sys.argv[1:])

 
抛开“ 包装”的代码,BBP(N)的功能是你真正想要的。你给它越大的N和给 getcontext().prec 设置越大的值,你就会使计算越精确。让我们看看一些代码结果:

这有许多数字位。你可以看出,我们并没有比以前更准确。所以我们需要前进到下一个公式,贝拉公式,希望能获得更好的精度。它看起来像这样: 

 我们将只改变我们的变换公式,其余的代码将保持不变。点击这里下载Python实现的贝拉公式。让我们看一看bellards(n):
 

def bellard(n):
  pi=Decimal(0)
  k=0
  while k < n:
    pi+=(Decimal(-1)**k/(1024**k))*( Decimal(256)/(10*k+1)+Decimal(1)/(10*k+9)-Decimal(64)/(10*k+3)-Decimal(32)/(4*k+1)-Decimal(4)/(10*k+5)-Decimal(4)/(10*k+7)-Decimal(1)/(4*k+3))
    k+=1
  pi=pi*1/(2**6)
  return pi

   哦,不,我们得到的是同样的精度。好吧,让我们试试第三个公式, Chudnovsky 算法,它看起来是这个样子: 

   再一次,让我们看一下这个计算公式(假设我们有一个阶乘公式)。 点击这里可下载用 python 实现的 Chudnovsky 公式。

下面是程序和输出结果:
 

def chudnovsky(n):
  pi=Decimal(0)
  k=0
  while k < n:
    pi+=(Decimal(-1)**k)*(Decimal(factorial(6*k))/((factorial(k)**3)*(factorial(3*k)))*(13591409+545140134*k)/(640320**(3*k)))
    k+=1
  pi=pi*Decimal(10005).sqrt()/4270934400
  pi=pi**(-1)
  return pi

    所以我们有了什么结论?花哨的算法不会使机器浮点世界达到更高标准。我真的很期待能有一个比我们用求和公式时所能得到的更好的精度。我猜那是过分的要求。如果你真的需要用PI,就只需使用math.pi变量了。然而,作为乐趣和测试你的计算机真的能有多快,你总是可以尝试第一个计算出Pi的百万位或者更多位是几。

以上是小编为您精心准备的的内容,在的博客、问答、公众号、人物、课程等栏目也有的相关内容,欢迎继续使用右上角搜索按钮进行搜索python
python实现四则运算、深入浅出python、深入python、深入python3、深入理解python,以便于您获取更多的相关知识。

时间: 2024-09-09 02:39:36

由Python运算π的值深入Python中科学计算的实现_python的相关文章

Python通过递归遍历出集合中所有元素的方法_python

本文实例讲述了Python通过递归遍历出集合中所有元素的方法.分享给大家供大家参考.具体实现方法如下: 复制代码 代码如下: '''''通过递归遍历出集合中的所有元素 Created on 2013-9-29 @author: L.Eric '''  def print_List(list_nums):      for each_item in list_nums :           if isinstance(each_item,list):              print_Lis

python字典键值对的添加和遍历方法_python

添加键值对 首先定义一个空字典 >>> dic={} 直接对字典中不存在的key进行赋值来添加 >>> dic['name']='zhangsan' >>> dic {'name': 'zhangsan'} 如果key或value都是变量也可以用这种方法 >>> key='age' >>> value=30 >>> dic[key]=value >>> dic {'age': 30

Python中用于计算对数的log()方法

 这篇文章主要介绍了Python中用于计算对数的log()方法,是Python入门基础中的必会的方法,需要的朋友可以参考下     log()方法返回x的自然对数,对于x>0. 语法 以下是log()方法的语法: ? 1 2 3 import math   math.log( x ) 注意:此函数是无法直接访问的,所以我们需要导入math模块,然后需要用math的静态对象来调用这个函数. 参数 x -- 这是一个数值表达式. 返回值 此方法返回x的自然对数,对于x>0. 例子 下面的例子显示了

详解Python中用于计算指数的exp()方法

  这篇文章主要介绍了详解Python中用于计算指数的exp()方法,是Python入门中必会的基本方法,需要的朋友可以参考下 exp()方法返回指数x: ex. 语法 以下是exp()方法的语法: ? 1 2 3 import math   math.exp( x ) 注意:此函数是无法直接访问的,所以我们需要导入math模块,然后需要用math的静态对象来调用这个函数. 参数 x -- 这是一个数值表达式 返回值 此方法返回指数x: ex. 例子 下面的例子显示了exp()方法的使用. ?

Python 和 Pandas 在 SQLite 数据库中的运用

本文讲的是Python 和 Pandas 在 SQLite 数据库中的运用, SQLite 是一个数据库引擎,使用它能方便地存储和处理关系型数据.它和 csv 格式很相似,SQLite 把数据存储在一个单独的文件中,它能方便地与其他人共享.大部分的编程语言和编译环境都对 SQLite 数据库提供了很好的支持.Python 也不例外,并且专门有一个访问 SQLite 数据库的程序库,叫做 sqlite3,自从 2.5 版本以来,它就已经被 Python 纳入标准库里.在这篇博文中,我们将学会如何使

浅析python 中__name__ = &#039;__main__&#039; 的作用_python

很多新手刚开始学习python的时候经常会看到python 中__name__ = \'__main__\' 这样的代码,可能很多新手一开始学习的时候都比较疑惑,python 中__name__ = '__main__' 的作用,到底干嘛的? 有句话经典的概括了这段代码的意义: "Make a script both importable and executable" 意思就是说让你写的脚本模块既可以导入到别的模块中用,另外该模块自己也可执行. 这句话,可能一开始听的还不是很懂.下面

Python中的Numpy入门教程_python

1.Numpy是什么 很简单,Numpy是Python的一个科学计算的库,提供了矩阵运算的功能,其一般与Scipy.matplotlib一起使用.其实,list已经提供了类似于矩阵的表示形式,不过numpy为我们提供了更多的函数.如果接触过matlab.scilab,那么numpy很好入手. 在以下的代码示例中,总是先导入了numpy: 复制代码 代码如下: >>> import numpy as np>>> print np.version.version1.6.2

python中元类用法实例_python

本文实例讲述了python中元类用法,分享给大家供大家参考.具体方法分析如下: 1.元类(metaclass)是用来创建类的类 2.type(object):返回一个对象的类型,与object.__class__的值相同,type(name,bases,dict):创建一个新的type类型,name就是新class的name,值存到__name__属性中,bases是tuple类型,值会存到__bases__中,dict的值存到__dict__中 复制代码 代码如下: class X: ... 

Python中list初始化方法示例_python

本文实例讲述了Python中list初始化方法.分享给大家供大家参考,具体如下: 1.基本方法. lst = [1, 2, 3, 4, 5] 2.初始化连续数字. >>> lst = [n for n in range(5, 10)] >>> print(lst) [5, 6, 7, 8, 9] 3.初始化n个相同值.(两种方式) >>> lst = ['x' for n in range(5)] >>> print(lst) ['x