沉降/引數識別

儘管努力使本篇關於絮凝懸浮液的引數識別儘可能自成一體，但對數值方法和懸浮液的建模的預備知識將會有用。特別是，牛頓-拉夫森方案用於解決非線性方程組的最佳化和有限體積法用於解決偏微分方程。

絮凝懸浮液的建模

絮凝懸浮液的分批沉降過程被建模為一個初值問題

${\frac {\partial u}{\partial t}}+{\frac {\partial f(u)}{\partial x}}={\frac {\partial A^{2}(u)}{\partial x^{2}}},\quad u(0,x)=u_{0}(),\quad x\in [0,L]$

其中$u$表示分散固相的體積分數。對於封閉性，對流通量函式由 Kynch 分批沉降函式給出，該函式具有 Richardson-Zaki 阻力函式

${\displaystyle f(u)=\left\{{\begin{matrix}&u_{\infty }u(1-u)^{C}&{\mbox{for }}0\leq u$

擴散通量由以下給出

$A(u)=\int _{0}^{u}a(s)ds,\quad a(u)=a_{0}(1-u)^{C}u^{k-1},$

它源於冪律的插入

$\sigma _{\mathrm {e} }(u)=\left\{{\begin{matrix}\sigma _{0}\left(\left({\frac {u}{u_{\mathrm {c} }}}\right)^{k}-1\right)&{\mbox{for }}u_{\mathrm {c} }<u\\0&{\mbox{for }}u\leq u_{\mathrm {c} }\end{matrix}}\right.$

到

$a(u)={\frac {f(u)\sigma _{\mathrm {e} }'(u)}{\Delta \varrho gu}},\quad a_{0}={\frac {u_{\infty }\sigma _{0}k}{\Delta \varrho gu_{\mathrm {c} }^{k}}}.$

在閉包中，常數$u_{\infty },C,a_{0},u_{\max },k$部分已知。

數值方案

直接問題的求解數值方案以守恆形式寫成內部點的行進公式（“內部方案”），如下所示：$u_{j}^{n+1}=u_{j}^{n}-{\frac {\Delta t}{\Delta x}}\left(f_{j+1/2}^{n+1}-f_{j-1/2}^{n+1}\right)+{\frac {\Delta t}{\Delta x^{2}}}\left(A(u_{j-1}^{n+1})-2A(u_{j}^{n+1})+A(u_{j+1}^{n+1})\right),\quad j=2,\dots ,J-1$

在邊界處（“邊界方案”）為

$u_{1}^{n+1}=u_{1}^{n}-{\frac {\Delta t}{\Delta x}}f_{j+1/2}^{n+1}+{\frac {\Delta t}{\Delta x^{2}}}\left(A(u_{j+1}^{n+1})-A(u_{j}^{n+1})\right),\quad u_{J}^{n+1}=u_{J}^{n}+{\frac {\Delta t}{\Delta x}}f_{J-1/2}^{n+1}-{\frac {\Delta t}{\Delta x^{2}}}\left(A(u_{J}^{n+1})-A(u_{J-1}^{n+1})\right)$

對於一階方案，數值通量函式變為

$f_{j+1/2}^{n}=f(u_{j}^{n},u_{j+1}^{n}).$ 如果通量函式只有一個最大值，表示為

$u_{\mathrm {m} }$，Engquist-Osher 數值通量函式可以表示為

$f^{\mathrm {EO} }(u,v)=\left\{{\begin{matrix}&f(u),&{\mbox{for }}u\leq u_{\mathrm {m} },v\leq u_{\mathrm {m} },\\&f(u)+f(v)-f(u_{\mathrm {m} }),&{\mbox{for }}u\leq u_{\mathrm {m} },v>u_{\mathrm {m} },\\&f(u_{\mathrm {m} }),&{\mbox{for }}u>u_{\mathrm {m} },v\leq u_{\mathrm {m} },\\&f(v),&{\mbox{for }}u>u_{\mathrm {m} },v>u_{\mathrm {m} }.\end{matrix}}\right.$

為了線性化，泰勒 公式

$f(u_{j}^{n+1},u_{j+1}^{n+1})=f(u_{j}^{n},u_{j+1}^{n})+{\frac {\partial f(u_{j}^{n},u_{j+1})}{\partial u_{j}^{n}}}(u_{j}^{n+1}-u_{j}^{n})+{\frac {\partial f(u_{j}^{n},u_{j+1})}{\partial u_{j+1}^{n}}}(u_{j+1}^{n+1}-u_{j+1}^{n}),\quad j=1,\dots ,J$

和

$A(u_{j}^{n+1})=A(u_{j}^{n})+{\frac {\partial A(u_{j}^{n})}{\partial u_{j}^{n}}}(u_{j}^{n+1}-u_{j}^{n}),\quad j=1,\dots ,J$

被插入。將時間演化步驟簡寫為

$\Delta u_{j}^{n}=u_{j}^{n+1}-u_{j}^{n},\quad j=1,\dots ,J,$

內部方案的線性化行進公式變為

$\Delta u_{j}^{n}=-{\frac {\Delta t}{\Delta x}}\left(f_{j+1/2}^{n}+{\mathcal {J}}_{j+1/2}^{-}\Delta u_{j}^{n}+{\mathcal {J}}_{j+1/2}^{+}\Delta u_{j+1}^{n}+f_{j-1/2}^{n}+{\mathcal {J}}_{j-1/2}^{-}\Delta u_{j-1}^{n}+{\mathcal {J}}_{j-1/2}^{+}\Delta u_{j}^{n}\right)$ $+{\frac {\Delta t}{\Delta x^{2}}}\left(A(u_{j-1})-2A(u_{j})+A(u_{j+1})+a(u_{j-1})\Delta u_{j-1}-2a(u_{j})\Delta u_{j}+a(u_{j+1})\Delta u_{j+1}\right),$

其中

$J_{j+1/2}^{+}={\frac {\partial f(u_{j}^{n},u_{j+1}^{n})}{\partial u_{j+1}^{n}}},\quad J_{j+1/2}^{-}={\frac {\partial f(u_{j}^{n},u_{j+1}^{n})}{\partial u_{j}^{n}}}.$

重新排列會得到一個塊三角線性系統

$-\left({\frac {\Delta t}{\Delta x}}({\mathcal {J}}_{j-1/2}^{-}+{\frac {\Delta t}{\Delta x^{2}}}a(u_{j-1}^{n})\right)\Delta u_{j-1}^{n}+\left(I+{\frac {\Delta t}{\Delta x}}\left({\mathcal {J}}_{j+1/2}^{-}-{\mathcal {J}}_{j-1/2}^{+}\right)+{\frac {2\Delta t}{\Delta x^{2}}}a(u_{j}^{n})\right)\Delta u_{j}^{n}$ $+\left({\frac {\Delta t}{\Delta x}}{\mathcal {J}}_{j+1/2}^{+}-{\frac {\Delta t}{\Delta x^{2}}}a(u_{j+1}^{n})\right)\Delta u_{j+1}^{n}$ $={\frac {\Delta t}{\Delta x}}\left(f(u_{j}^{n},u_{j+1}^{n})-f(u_{j-1}^{n},u_{j}^{n})\right)-{\frac {\Delta t}{\Delta x^{2}}}\left(A(u_{j-1}^{n})-2A(u_{j}^{n})+A(u_{j+1}^{n})\right),j=2,\dots ,J-2$

這是一種形式

$m_{j,j-1}\Delta u_{j-1}^{n}+m_{j,j}\Delta u_{j}^{n}+m_{j,j+1}\Delta u_{j+1}=b_{j},$

或者，用更緊湊的符號，

$M\Delta u^{n}=b,\quad {\begin{bmatrix}m_{11}&m_{12}&0&\cdots &\cdots &\cdots &\cdots &0\\m_{21}&m_{22}&m_{23}&0&\cdots &\cdots &\cdots &0\\\vdots &&&&&&&\vdots \\\vdots &&&&&&\ddots &\vdots \\0&\cdots &\cdots &\cdots &\cdots &0&m_{J,J-1}&m_{J,J}\end{bmatrix}},\quad b={\begin{bmatrix}b_{1}\\\vdots \\b_{J}\end{bmatrix}},\quad \Delta u^{n}={\begin{bmatrix}\Delta u_{1}^{n}\\\vdots \\\Delta u_{J}^{n}\end{bmatrix}}.$

引數識別作為最佳化問題

透過最佳化方法進行引數識別，目標是在引數空間上最小化成本函式

$\min _{e}q(e),\quad q(e)=\|h(e)-H\|^{2}$,

其中 h(e) 表示透過模擬計算的介面，H 表示測量的介面。不失一般性，我們考慮一個由兩個引數組成的引數集 e=(e_1, e_2)。最佳化可以透過使用牛頓法迭代地進行，如下所示

$Q(e_{k})\Delta e_{k}=\nabla _{e}q(e_{k}),\quad e\subset \{C,v_{\infty },a_{0},u_{\mathrm {c} },k\},$

其中

$Q={\begin{bmatrix}{\frac {\partial ^{2}q}{\partial e_{1}^{2}}}&{\frac {\partial ^{2}q}{\partial e_{1}e_{2}}}\\{\frac {\partial ^{2}q}{\partial e_{2}e_{1}}}&{\frac {\partial ^{2}q}{\partial e_{2}^{2}}}\end{bmatrix}},$

是 q 的海森矩陣。牛頓法是基於泰勒展開式

$0\approx Q(e^{*})=Q(e)+\nabla _{e}Q(e)\Delta e+{\frac {1}{2}}\Delta e^{\mathrm {T} }Q\Delta e,\quad \Delta e=e^{*}-e,$

其中 $e^{*}$ 是最佳的引數選擇。

外部連結

• http://de.wikipedia.org/wiki/Wikipedia:TeX 維基百科：TeX
• http://en.wikipedia.org/wiki/Matlab