# Demo entry 6356246

important sampling

Submitted by anonymous on Apr 17, 2017 at 17:55
Language: Python 3. Code size: 2.2 kB.

```# -*- coding: utf-8 -*-
"""
Created on Wed Apr 12 14:29:10 2017

@author: think
"""

import numpy as np
import matplotlib.pyplot as plt

AA=0.0  #左边界
BB=20.0 #右边界
width=BB-AA  #区间长度
np.random.seed(456) #随机数种子

#所要求解的被积函数表达式
def ff(x):
return x**4*np.exp(-x)

#权重函数1
def weight1(x):
return 1.0/width

#权重函数2
def weight2(x):
return np.exp(-x)

#权重函数3
def weight3(x):
return x**3*np.exp(-x)/6.0

#权重函数4
def weight4(x):
return x**4*np.exp(-x)/24.0

maxSave=[1/width,1.0,0.3,0.3] #每个权重函数的最大值
#每个权重函数的引用
weightSave=[weight1,weight2,weight3,weight4]

"""

"""
cal_num=25 #执行的总的计算次数
results1=[] #储存权重函数1情形的结果
results3=[] #储存权重函数3的情形的结果
results=[[],[]]
chooseIndex=[1,3] #即选择1和3两种情形
nt=2000  #每次计算时的样本容量

"""

"""
dd=0
while dd<np.size(chooseIndex):
kk=chooseIndex[dd]
mm=0
while mm<cal_num:
"""
mylist是用于储存所有的可以作为样本元素
的delta的值，无需关心它们存入的先后顺序
"""
mylist=[]
i=0
#产生符合指定分布的样本
while i<nt:
delta=AA+np.random.random()*width
fvalue=weightSave[kk-1](delta)/maxSave[kk-1]
randJudge=np.random.random()
if randJudge<fvalue:
mylist.append(delta)
i+=1

#按重要性取样计算积分值
jj=0.0
for ii in range(np.size(mylist)):
x=mylist[ii]
jj+=ff(x)/weightSave[kk-1](x)
jj=jj/np.size(mylist)
results[dd].append(jj) #储存结果
mm+=1
dd+=1

results1=results[0]
results3=results[1]

mean1=round(np.mean(results1),2)
std1=round(np.std(results1),2)
sterr1=round(std1/np.sqrt(cal_num),2)

mean3=round(np.mean(results3),2)
std3=round(np.std(results3),2)
sterr3=round(std3/np.sqrt(cal_num),2)

print('权重1对应的平均值',mean1,'权重3对应的平均值',mean3)
print('权重1对应的标准差',std1,'权重3对应的标准差',std3)
print('权重1对应的标准误',sterr1,'权重3对应的标准误',sterr3)

plt.plot(results1,'r',label='weight1')
plt.plot(results3,'b',label='weight3')
plt.ylim(20,28)
plt.xlabel('X')
plt.ylabel('Y')
plt.legend(loc='upper right')
plt.show()
```

This snippet took 0.00 seconds to highlight.

Back to the Entry List or Home.