有限元方法入門編程

我是張運霄，歡迎來到學習的最優體驗

解微分方程有很多數值方法，有限差分，有限元，有限體積，邊界元等，它們本質上是用不同的方式取離散化微分方程。今天說一下有限元方法入門編程。

  問題描述  

我們有如下問題，要求解u(x)。理論部分如下，選定基函數，基函數屬於勒貝格空間(見附錄)，我們對目標微分方程兩邊乘以測試函數後分部積分，得到弱格式，最終得到一個線性系統，求解線性方程即可得到微分方程的近似解。

  編程實現（只核心代碼） 

1 剛度矩陣A計算

for n in range(N):    E = Tb[:,n]    for alpha in range(Nlb):        for beta in range(Nlb):            ff = lambda x: c(x)*(-1)**(alpha+beta)/h**2            A[E[beta],E[alpha]] += quad(ff,P[E[0]],P[E[1]])[0]
2 力向量b計算
for n in range(N):    E = Tb[:,n]    for beta in range(Nlb):        ff = lambda l: f(l)*psi(l,E,beta)        b[E[beta]] += quad(ff,P[E[0]],P[E[1]])[0]
3 邊界條件處理
for k in range(Nbn):    if boundary[0][k] == 'Dirichlet':        i = boundary[1][k]        A[i,:]=0        A[i,i]=1        b[i]=g(Pb[i])
4 求解與誤差計算
u = solve(A, b)u_exact = guu = u_exact(Pb)error = max(abs(u - uu))print(error)
5 作圖
  附錄：有限元理論中涉及到的幾個概念  
支集
一個函數的支集是它定義域裡使得函數不為零的範圍構成的子集
例如：delta(x)的支集為{0}，因為它只在0處等於1
    
緊支集
有邊界的支集
[1,3]為緊支集，[1,+inf)不是
函數可導性分類
C^0 所有連續函數
C^1 一階可導且導數連續的集合
C^k 所有k階可導且導數連續的函數集合
C^inf 無限可導
弱導數
如果下式成立(int為積分)
int{a,b}{wv} = -int{a,b}{uv'}
則w是u的弱導數
勒貝格積分 Lebesgue integral
黎曼積分把積分區域縱向離散
勒貝格積分把積分區域橫向離散
勒貝格空間 Lebesgue spaces，L^p spaces
有限元空間與p-範數構成的空間（p-範數不為無窮）

p-範數
def norm_f(f,a,b,p):    """    ||x||_p = int{a,b}{|f|**p}**(1/p)    """        if p == "inf":        fp = lambda x: np.abs(f(x))        return fmin(fp)    fp = lambda x: np.abs(f(x))**p    return quad(fp,a,b)[0]**(1/p)

半賦范 seminorm
回顧範數norm的性質如下
commutes with scaling
obeys triangle inequality
zero only at the origin
半賦范少了範數的第三條，也就是不在原點的向量也可能範數為0。
像這樣：x=[1,0,0], 但  ||x|| = 0

賦范向量空間 normed vector spaces
向量V和範數||·||構成的空間(V,||V||)

半賦范向量空間 seminormed vector spaces
向量V和範數p-norm構成的空間(V,p)

內積空間
向量V和內積<,>構成的空間(V,<,>)
內積空間是歐幾裡得空間(三維)的擴展(到無窮維)
歐幾裡得空間中的內積擴展為點乘

度量空間 metric space
一個集合M及其度量d(一般是個函數)構成的空間(M,d)

巴拿赫空間 Banach space
巴拿赫空間是完備賦范向量空間 complete normed vector space

索博列夫空間 Sobolev space
向量V和L^p範數構成的空間
W(k,p),k表徵導數性質，p表徵範數性質

希爾伯特空間 Hilbert space
當p=2時
H^k = W(k,2)

  參考  
對有限元編程感興趣的可以去b站看何曉明的有限元編程短課，每周末直播。
https://en.wikipedia.org/wiki/Function_space
長按關注學習的最優體驗收穫更多