python实现蒙特卡罗模拟法的实践

 更新时间:2022年03月04日 09:43:38   作者:洋洋菜鸟  
 蒙特卡洛就是产生随机变量,带入模型算的结果,寻优方面,本文主要介绍了python 蒙特卡罗模拟法实现,文中通过示例代码介绍的非常详细,具有一定的参考价值,感兴趣的小伙伴们可以参考一下

1.简介

        蒙特卡洛又称随机抽样或统计试验,就是产生随机变量,带入模型算的结果,寻优方面,只要模拟次数够多,最终是可以找到最优解或接近最优的解。

         通常蒙特卡罗方法可以粗略地分成两类:一类是所求解的问题本身具有内在的随机性,借助计算机的运算能力可以直接模拟这种随机的过程。例如在核物理研究中,分析中子在反应堆中的传输过程。中子与原子核作用受到量子力学规律的制约,人们只能知道它们相互作用发生的概率,却无法准确获得中子与原子核作用时的位置以及裂变产生的新中子的行进速率和方向。科学家依据其概率进行随机抽样得到裂变位置、速度和方向,这样模拟大量中子的行为后,经过统计就能获得中子传输的范围,作为反应堆设计的依据。
另一种类型是所求解问题可以转化为某种随机分布的特征数,比如随机事件出现的概率,或者随机变量的期望值。通过随机抽样的方法,以随机事件出现的频率估计其概率,或者以抽样的数字特征估算随机变量的数字

2.实例分析

2.1 模拟求近似圆周率

        绘制单位圆和外接正方形,正方形ABCD的面积为:2*2=4,圆的面积为:S=Π*1*1=Π,现在模拟产生在正方形ABCD中均匀分布的点n个,如果这n个点中有m个点在该圆内,则圆的面积与正方形ABCD的面积之比可近似为m/n

程序如下:

#模拟求近似圆周率
import random
import numpy as np
import matplotlib.pyplot as plt
#解决图标题中文乱码问题
import matplotlib as mpl
mpl.rcParams['font.sans-serif'] = ['SimHei']  # 指定默认字体
mpl.rcParams['axes.unicode_minus'] = False  # 解决保存图像是负号'-'显示为方块的问题
#进入正题
r=random.random()
num=range(0,200000,10)
mypi=np.ones((1,len(num)))
for j in range(len(num)):
    # 投点次数
    n = 10000
    # 圆的半径、圆心
    r = 1.0
    a,b = (0.,0.)
    # 正方形区域
    x_min, x_max = a-r, a+r
    y_min, y_max = b-r, b+r
    # 在正方形区域内随机投点
    x = np.random.uniform(x_min, x_max, n) #均匀分布
    y = np.random.uniform(y_min, y_max, n)
    #计算点到圆心的距离
    d = np.sqrt((x-a)**2 + (y-b)**2)
    #统计落在圆内点的数目
    res = sum(np.where(d < r, 1, 0))
    #计算pi的近似值(Monte Carlo:用统计值去近似真实值)
    mypi[0,j] = 4 * res / n
plt.plot(range(1,len(mypi[0])+1),mypi[0],'.-') 
plt.grid(ls=":",c='b',)#打开坐标网格
plt.axhline(y=np.pi,ls=":",c="yellow")#添加水平直线
# plt.axvline(x=4,ls="-",c="green")#添加垂直直线
plt.legend(['模拟', '实际'], loc='upper right', scatterpoints=1)
plt.title("近似圆周率")
plt.show()

返回:

2.2 估算定积分

程序如下:

#估算定积分
import random
import numpy as np
import matplotlib.pyplot as plt
#解决图标题中文乱码问题
import matplotlib as mpl
mpl.rcParams['font.sans-serif'] = ['SimHei']  # 指定默认字体
mpl.rcParams['axes.unicode_minus'] = False  # 解决保存图像是负号'-'显示为方块的问题
#进入正题
r=random.random()
num=range(1,10**6,500)
s = np.ones((1,len(num)))
for j in range(len(num)):
    n = 30000
    #矩形边界区域
    x_min, x_max = 0.0, 1.0
    y_min, y_max = 0.0, 1.0
    #在矩形区域内随机投点x
    x = np.random.uniform(x_min, x_max, n)#均匀分布
    y = np.random.uniform(y_min, y_max, n)
    #统计落在函数y=x^2下方的点
    res = sum(np.where(y < x**2, 1 ,0))
    #计算定积分的近似值
    s[0,j] = res / n
plt.plot(range(1,len(s[0])+1),s[0],'.-') 
plt.grid(ls=":",c='b',)#打开坐标网格
plt.axhline(y=1/3,ls=":",c="red")#添加水平直线
# plt.axvline(x=4,ls="-",c="green")#添加垂直直线
plt.legend(['模拟', '实际1/3'], loc='upper right', scatterpoints=1)
plt.title("估算定积分")
plt.show()

返回:

2.3 求解整数规划

要解的方程为:

条件如下:

程序如下:

# 求解整数规划
import random
import numpy as np
import time
import matplotlib.pyplot as plt
#解决图标题中文乱码问题
import matplotlib as mpl
mpl.rcParams['font.sans-serif'] = ['SimHei']  # 指定默认字体
mpl.rcParams['axes.unicode_minus'] = False  # 解决保存图像是负号'-'显示为方块的问题
#进入正题
time_start=time.time() #计时开始
p0=0
for i in range(10**7):
    x=np.random.randint(0,99,(1,3))
    f=2*x[0,0]+3*x[0,0]**2+3*x[0,1]+x[0,1]**2+x[0,2]
    g=[
        x[0,0]+2*x[0,0]**2+x[0,1]+2*x[0,1]**2+x[0,2],
        x[0,0]+x[0,0]**2+x[0,1]+x[0,1]**2-x[0,2],
        2*x[0,0]+x[0,0]**2+2*x[0,1]+x[0,2],
        x[0,0]+2*x[0,1]**2
    ]
    if g[0]<=100 and g[1]<=500 and g[2]<=400 and g[3]>=10:
        if p0<f:
            x0=x
            p0=f
print('最优解:',x0)
print('最优值:',p0)
time_end=time.time() #计时结束
print('用时:',time_end-time_start)

返回:

 到此这篇关于python实现蒙特卡罗模拟法的实践的文章就介绍到这了,更多相关python 蒙特卡罗模拟法内容请搜索脚本之家以前的文章或继续浏览下面的相关文章希望大家以后多多支持脚本之家!

相关文章

  • python调用支付宝支付接口流程

    python调用支付宝支付接口流程

    这篇文章主要介绍了python调用支付宝支付接口流程,本文给大家介绍的非常详细,具有一定的参考借鉴价值,需要的朋友可以参考下
    2019-08-08
  • 使用python解析MDX词典数据并保存为Excel文件

    使用python解析MDX词典数据并保存为Excel文件

    MDX(Mobile Dictionary eXchange)是一种常见的词典文件格式,通常用于在移动设备和电脑之间共享辞典数据,本文深入探讨了从MDX词典数据提取、处理到最终保存为Excel文件的全过程,需要的朋友可以参考下
    2023-12-12
  • 使用Python编写vim插件的简单示例

    使用Python编写vim插件的简单示例

    这篇文章主要介绍了使用Python编写vim插件的简单教程,文中举了一个获取reddit首页信息并显示在缓冲区中的例子,需要的朋友可以参考下
    2015-04-04
  • Python实现AES加密,解密的两种方法

    Python实现AES加密,解密的两种方法

    这篇文章主要介绍了Python实现AES加密,解密的两种方法,帮助大家更好的使用python加解密文件,感兴趣的朋友可以了解下
    2020-10-10
  • Python区块链客户端类开发教程

    Python区块链客户端类开发教程

    这篇文章主要为大家介绍了Python区块链客户端类开发教程,有需要的朋友可以借鉴参考下,希望能够有所帮助,祝大家多多进步,早日升职加薪
    2022-05-05
  • 在Pycharm中自动添加时间日期作者等信息的方法

    在Pycharm中自动添加时间日期作者等信息的方法

    今天小编就为大家分享一篇在Pycharm中自动添加时间日期作者等信息的方法,具有很好的参考价值,希望对大家有所帮助。一起跟随小编过来看看吧
    2019-01-01
  • python本地文件服务器实例教程

    python本地文件服务器实例教程

    这篇文章主要给大家介绍了关于python本地文件服务器的相关资料,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的朋友们下面随着小编来一起学习学习吧
    2021-05-05
  • 简单了解Python3 bytes和str类型的区别和联系

    简单了解Python3 bytes和str类型的区别和联系

    这篇文章主要介绍了简单了解Python3 bytes和str类型的区别和联系,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的朋友可以参考下
    2019-12-12
  • PIP和conda 更换国内安装源的方法步骤

    PIP和conda 更换国内安装源的方法步骤

    这篇文章主要介绍了PIP和conda 更换国内安装源的方法步骤,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的朋友们下面随着小编来一起学习学习吧
    2020-09-09
  • Python与Matlab混合编程的实现案例

    Python与Matlab混合编程的实现案例

    本文主要介绍了Python与Matlab混合编程的实现案例,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的朋友们下面随着小编来一起学习学习吧
    2023-01-01

最新评论