mod_ppm Module

The one-dimensional implementation of the PPM （piecewise parabolic method, 分段抛物线法）

• 假设1：等距网格
• 假设2：第一个网格和最后一个网格为区域边界网格，为外部约束。

Step 1: 计算2~n-1个网格的平均坡度(slopC)

$$\delta c_{i} = 0.5(c_{i+1} - c_{i-1})$$

为了让浓度不连续的地方模拟得更好， 同时，确保网格边界处浓度在左右格子平均浓度之间， 对$\delta c_{i}$进行约束。 $$\delta_m c_{i} = \begin{cases} min(|\delta c_{i}|, 2|c_{i+1} - c_{i}|, 2|c_{i} - c_{i-1}|)*sign(\delta c_{i}) &\quad {\text{if } (c_{i+1} - c_{i})(c_{i} - c_{i-1}) \gt 0 }\\ \text{0,} &\quad \text{else} \\ \end{cases} \tag {2.9} \label {2.9}$$

Step 2: 每个网格边界处的浓度 $$c_{x_{i+1/2}} = c_i + 0.5(c_{i+1}-c_i) + \frac{1}{6}(\delta_m c_{i} - \delta_m c_{i+1}) \tag {2.8} \label {2.8}$$ Step 3: 计算抛物线参数，并对边界处浓度进行约束 $$c_{6,i} = 6(c_i- \frac{1}{2}(c_{x_{i+1/2}} + c_{x_{i-1/2}}) ) \tag {2.5} \label {2.5}$$ $$\begin{align*} c_{x_{i-1/2}} &= c_{x_{i+1/2}} = c_{i} \quad \text{if } (c_{x_{i+1/2}}-c_{i}) (c_{i}-c_{x_{i-1/2}}) \le 0 \\ c_{x_{i-1/2}} &= 3c_{i} - 2c_{x_{i+1/2}} \quad \text{if } \Delta c_i c_{6,i} > (\Delta c_i)^2 \\ c_{x_{i+1/2}} &= 3c_{i} - 2c_{x_{i-1/2}} \quad \text{if } - \Delta c_i c_{6,i} > (\Delta c_i)^2 \tag {2.10} \label {2.10} \end{align*}$$

Step 4: 计算边界处传输的平均浓度，及其通量 $$\begin{align*} c_{i+1/2, L}(y) &= c_{x_{i+1/2}} - \frac{y}{2\Delta x} (c_{x_{i+1/2}} - c_{x_{i-1/2}} -(1-\frac{2y}{3\Delta x})c_{6,i}) \\ c_{i+1/2, R}(y) &= c_{x_{i+1-1/2}} - \frac{y}{2\Delta x} (c_{x_{i+1+1/2}} - c_{x_{i+1-1/2}} -(1-\frac{2y}{3\Delta x})c_{6,i+1}) \tag {2.12} \label {2.12} \end{align*}$$ 其中$y=|u\Delta t|$ $$f_{i+1/2} = \begin{cases} u_{x_{i+1/2}} \times c_{i+1/2, L} &\quad {\text{if } u_{x_{i+1/2}} \ge 0 }\\ u_{x_{i+1/2}} \times c_{i+1/2, R} &\quad {\text{if } u_{x_{i+1/2}} \le 0 } \\ \end{cases} \tag {2.14} \label {2.14}$$ Step 5: 用前向欧拉法进行时间积分 $$c_i^{n+1} = c_i^{n} + \frac{\Delta t}{\Delta x}(f_{i-1/2}-f_{i+1/2}) \tag {2.13} \label {2.13}$$ 网格设计

     |Boundary|<-----------------Domain----------------->|Boundary|
     |  c(1)  |   c(2)  |-----------------------| c(n-1) |  c(n)  |     
     |       u(1)      u(2) ->                         u(n-1)    u(n)
     |        | slopC(2)|
     |        |deltaC(2)|
     |        |  c6(2)  |
  leftC(1) leftC(2) 
           rightC(1) rightC(2) 
  egdeC(0) egdeC(1)  egdeC(2)
           transC(1) transC(2)
           flux(1)    flux(2) ->                     flux(n-1)