(本文經Prof. Lorena A. Barba授權轉載翻譯,原文鏈結請參照底下參考資料)
CFD Python:求解Navier-Stokes方程的12個步驟(1) CFD Python:求解Navier-Stokes方程的12個步驟(2) CFD Python:求解Navier-Stokes方程的12個步驟(3) CFD Python:求解Navier-Stokes方程的12個步驟(4) CFD Python:求解Navier-Stokes方程的12個步驟(5) CFD Python:求解Navier-Stokes方程的12個步驟(6) CFD Python:求解Navier-Stokes方程的12個步驟(7) CFD Python:求解Navier-Stokes方程的12個步驟(8)
步驟2:一維非線性對流方程
現在要用步驟1的方法來實踐非線性對流方程,一維對流方程表示如下:
\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = 0
相較步驟1的方程式第二項乘以常數c,現在第二項乘以u。因此,方程式第二項為非線性。使用如同步驟1的離散方法,即時間項的前向差分與空間項的後向差分來離散。離散方程如下:
\frac{u_i^{n+1}-u_i^n}{\Delta t} + u_i^n \frac{u_i^n-u_{i-1}^n}{\Delta x} = 0
為了求解未知項u_i^{n+1},可將方程式改寫如下:
u_i^{n+1} = u_i^n - u_i^n \frac{\Delta t}{\Delta x} (u_i^n - u_{i-1}^n)
如同之前作法,Python方程式需先引入必要的函式庫,接著宣告一些空間項與時間項的變數(你應該更改參數值來看看會發生什麼)。再來給定初始條件u_0來初始化存放待求變數的陣列。初始條件為u = 2\ at \ 0.5 \leq x \leq 1且(0, 2)其他任何地方u = 1(此為 hat function)。
import numpy # we're importing numpy
from matplotlib import pyplot # and our 2D plotting library
%matplotlib inline
nx = 41
dx = 2 / (nx - 1)
nt = 20 #nt is the number of timesteps we want to calculate
dt = .025 #dt is the amount of time each timestep covers (delta t)
u = numpy.ones(nx) #as before, we initialize u with every value equal to 1.
u[int(.5 / dx) : int(1 / dx + 1)] = 2 #then set u = 2 between 0.5 and 1 as per our I.C.s
un = numpy.ones(nx) #initialize our placeholder array un, to hold the time-stepped solution
下方程式碼還沒有完成,只是單純從步驟1複製程式碼來執行時間步階的迭代,你能編輯這段程式碼來求解非線性對流方程嗎?
for n in range(nt): #iterate through time
un = u.copy() ##copy the existing values of u into un
for i in range(1, nx): ##now we'll iterate through the u array
###This is the line from Step 1, copied exactly. Edit it for our new equation.
###then uncomment it and run the cell to evaluate Step 2
###u[i] = un[i] - c * dt / dx * (un[i] - un[i-1])
pyplot.plot(numpy.linspace(0, 2, nx), u) ##Plot the results
你在非線性對流方程下從帽函式(hat function)的演變觀察到什麼呢?當你更改參數值且重新執行後又發生了什麼呢?
參考資料
- 12 steps to Navier–Stokes: Step 2: Nonlinear Convection,Lorena A. Barba