蒙特卡罗模拟法 —— python
生活随笔
收集整理的這篇文章主要介紹了
蒙特卡罗模拟法 —— python
小編覺得挺不錯的,現在分享給大家,幫大家做個參考.
目錄
1.簡介
2.實例分析
2.1 模擬求近似圓周率
2.2?估算定積分
2.3?求解整數規劃
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.0a,b = (0.,0.)# 正方形區域x_min, x_max = a-r, a+ry_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.0y_min, y_max = 0.0, 1.0#在矩形區域內隨機投點xx = 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=xp0=f print('最優解:',x0) print('最優值:',p0) time_end=time.time() #計時結束 print('用時:',time_end-time_start)返回:
總結
以上是生活随笔為你收集整理的蒙特卡罗模拟法 —— python的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: oracle11 处理器数,11G AW
- 下一篇: CRITIC法之python