通过推广矩形的方法整合3维
integration 3 dimensions by generalizing the method of rectangle
def RectanglesPointMilieu(f,a,b,n):
interval = 1.* (b-a) / n
sumfct = 0
for i in np.arange(a,b,interval):
sumfct += f(i + interval/2.)
return interval * sumfct
如何将其更改为 3 维计算?
def rec3(f,X1,X2,X3,a1,a2,a3,b1,b2,b3,N1,N2,N3):
interval1 = 1.* (b1-a1) / N1
interval2 = 1.* (b2-a2) / N2
interval3 = 1.* (b3-a3) / N3
Sum1 = 0
Sum2 = 0
Sum3 = 0
for i in np.arange(a1,b1,interval1):
Sum1 += f((X1[i]+ interval1))
for i in np.arange(a2,b2,interval2):
Sum2 += f((X2[i]+ interval2))
for i in np.arange(a3,b3,interval3):
Sum3 += f((X3[i]+ interval3))
return interval * float(Sum)
我做了类似的事情,但我是如此迷茫和困惑...我不知道如何继续...[我是新手]
def rec3(f,X1,X2,X3,a1,a2,a3,b1,b2,b3,N1,N2,N3):
interval1 = 1.* (b1-a1) / N1
interval2 = 1.* (b2-a2) / N2
interval3 = 1.* (b3-a3) / N3
Sum = 0
for i in np.arange(a1,b1,interval1):
for j in np.arange(a2,b2,interval2):
for k in np.arange(a3,b3,interval3):
Sum += f((X1[i]+ interval1),(X2[j]+ interval2),(X3[k]+ interval3))
return interval1 * interval2 * interval3 * float(Sum)
与
类型错误:'numpy.ndarray' 对象不可调用
#def TEST_Q2():
# Créer des tableaux
N = 1E5
X1 = rd.normal(0,1,N)
X2 = rd.normal(0,1,N)
X3 = rd.normal(0,1,N)
XX1 = X1 * np.sqrt(2)
XX2 = X2 * np.sqrt(2)
XX3 = X3 * 2.
def fct(a,b,c): # [x1**2 * x2**4 * exp(-x1**2 -x2**2 -2x3**2)]
return ((a**2. * b**4.)/32.) * np.exp(-a**2./2.-b**2./2.-c**2./2.)
# Calculer l’intégrale
F = fct(XX1,XX2,XX3)
print rec3(F,XX1,XX2,XX3,0,0,0,1,1,1,N,N,N)
# Rlt = 11.8416033988
我用这个作为我的测试
这叫做multidimensional numerical integration。为了在 n
维函数上执行此操作,您将需要使用 n
嵌套循环。在您 n=3
的情况下,天真的方法是简单地做类似的事情:
def rec3(f,a1,a2,a3,b1,b2,b3,N1,N2,N3):
interval1 = 1.* (b1-a1) / N1
interval2 = 1.* (b2-a2) / N2
interval3 = 1.* (b3-a3) / N3
Sum = 0
for x in np.arange(a1,b1,interval1):
for y in np.arange(a2,b2,interval2):
for z in np.arange(a3,b3,interval3):
Sum += f(x+interval1/2., y+interval2/2., z+interval3/2.)
return abs(interval1*interval2*interval3)*Sum
我不是很了解python,所以里面可能有语法错误,请谨慎使用。
还有其他方法可能会产生更好的结果。您可以尝试我在上面建议的搜索中的一些链接,看看是否有更适合您特定需求的链接。
更新:这是它的用法:
N = 1E3
def fct(a,b,c):
return ((a**2. * b**4.)/32.) * np.exp(-a**2./2.-b**2./2.-c**2./2.)
print rec3(fct,0,0,0,1,1,1,N,N,N)
def RectanglesPointMilieu(f,a,b,n):
interval = 1.* (b-a) / n
sumfct = 0
for i in np.arange(a,b,interval):
sumfct += f(i + interval/2.)
return interval * sumfct
如何将其更改为 3 维计算?
def rec3(f,X1,X2,X3,a1,a2,a3,b1,b2,b3,N1,N2,N3):
interval1 = 1.* (b1-a1) / N1
interval2 = 1.* (b2-a2) / N2
interval3 = 1.* (b3-a3) / N3
Sum1 = 0
Sum2 = 0
Sum3 = 0
for i in np.arange(a1,b1,interval1):
Sum1 += f((X1[i]+ interval1))
for i in np.arange(a2,b2,interval2):
Sum2 += f((X2[i]+ interval2))
for i in np.arange(a3,b3,interval3):
Sum3 += f((X3[i]+ interval3))
return interval * float(Sum)
我做了类似的事情,但我是如此迷茫和困惑...我不知道如何继续...[我是新手]
def rec3(f,X1,X2,X3,a1,a2,a3,b1,b2,b3,N1,N2,N3):
interval1 = 1.* (b1-a1) / N1
interval2 = 1.* (b2-a2) / N2
interval3 = 1.* (b3-a3) / N3
Sum = 0
for i in np.arange(a1,b1,interval1):
for j in np.arange(a2,b2,interval2):
for k in np.arange(a3,b3,interval3):
Sum += f((X1[i]+ interval1),(X2[j]+ interval2),(X3[k]+ interval3))
return interval1 * interval2 * interval3 * float(Sum)
与 类型错误:'numpy.ndarray' 对象不可调用
#def TEST_Q2():
# Créer des tableaux
N = 1E5
X1 = rd.normal(0,1,N)
X2 = rd.normal(0,1,N)
X3 = rd.normal(0,1,N)
XX1 = X1 * np.sqrt(2)
XX2 = X2 * np.sqrt(2)
XX3 = X3 * 2.
def fct(a,b,c): # [x1**2 * x2**4 * exp(-x1**2 -x2**2 -2x3**2)]
return ((a**2. * b**4.)/32.) * np.exp(-a**2./2.-b**2./2.-c**2./2.)
# Calculer l’intégrale
F = fct(XX1,XX2,XX3)
print rec3(F,XX1,XX2,XX3,0,0,0,1,1,1,N,N,N)
# Rlt = 11.8416033988
我用这个作为我的测试
这叫做multidimensional numerical integration。为了在 n
维函数上执行此操作,您将需要使用 n
嵌套循环。在您 n=3
的情况下,天真的方法是简单地做类似的事情:
def rec3(f,a1,a2,a3,b1,b2,b3,N1,N2,N3):
interval1 = 1.* (b1-a1) / N1
interval2 = 1.* (b2-a2) / N2
interval3 = 1.* (b3-a3) / N3
Sum = 0
for x in np.arange(a1,b1,interval1):
for y in np.arange(a2,b2,interval2):
for z in np.arange(a3,b3,interval3):
Sum += f(x+interval1/2., y+interval2/2., z+interval3/2.)
return abs(interval1*interval2*interval3)*Sum
我不是很了解python,所以里面可能有语法错误,请谨慎使用。
还有其他方法可能会产生更好的结果。您可以尝试我在上面建议的搜索中的一些链接,看看是否有更适合您特定需求的链接。
更新:这是它的用法:
N = 1E3
def fct(a,b,c):
return ((a**2. * b**4.)/32.) * np.exp(-a**2./2.-b**2./2.-c**2./2.)
print rec3(fct,0,0,0,1,1,1,N,N,N)