复合高斯积分(节点数小于等于3的版本Python实现)
生活随笔
收集整理的這篇文章主要介紹了
复合高斯积分(节点数小于等于3的版本Python实现)
小編覺得挺不錯的,現在分享給大家,幫大家做個參考.
被積函數
y=11+x2y=11+x2
x∈[?1,1]x∈[?1,1]
算法描述
缺點: 在下面代碼中,沒有實現點數大于3的情況(等于3的話可以使用)
因為計算對應的系數值的時候,這里就不是直接用python的sympy似乎就搞不定了。
簡述: 就是通過勒讓德(正交多項式)的零點。來做積分的考量點。之后再結合高斯積分的特點,得到對應的代數精度下方程,解出對應區間內的系數值。
當然也可以直接做伸縮變換。之后就需要乘以一個系數。(在積分的微分變量求導,導出的一個系數)
這個思想類似于之前做的在契比雪夫多項式零點的拉格朗日插值。
這些思想告訴我們,其實在做逼近之類的事情的時候,我們如果可以發現在某些特別重要的特征,重點將這些特征逼近,就可以大大的提高效率。類似于信息論里的不定長編碼,通過概率來導出最短信源編碼長一樣,來使得整體信源符號熵最小的過程類似。
都是強調對重要特征的重點描述。
代碼
import numpy as np from sympy import *def t(begin, end):xs = ps * (end - begin) / 2 + (begin + end) / 2temp = 0for i in range(point):temp += al[i] * Y.subs(x, xs[i])temp = temp.subs(ans)return temp.evalf() * (end - begin) / 2def P():if point == 0:return 1elif point == 1:return xp0 = 1p1 = xfor i in range(point - 1):temp = ((2 * i + 3) * x * p1 - (i + 1) * p0) / (i + 2)p0 = p1p1 = tempreturn p1def loss(begin=-1, end=1):T = sum([t(xl[i], xl[i + 1]) for i in range(n)])I = integrate(Y, (x, begin, end))print('%.18f' % (I - T).evalf())if __name__ == '__main__':x = symbols('x')point = 2n = 10xl = np.linspace(-1, 1, n + 1)Y = 1 / (1 + x ** 2)ps = np.array(solve(P()))funcs = []al = symbols('a0:%d' % point)for i in range(point):if i != 0:temp = integrate((x ** i), (x, -1, 1))for j in range(point):temp -= al[j] * (ps[j] ** i)else:temp = integrate(1, (x, -1, 1))for j in range(point):temp -= al[j]funcs.append(temp)ans = solve(funcs)loss()總結
以上是生活随笔為你收集整理的复合高斯积分(节点数小于等于3的版本Python实现)的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 复合五点高斯公式计算(Python实现)
- 下一篇: 多线程读取矩阵文件+多线程矩阵乘法(C+