计算科学 - 使用 FEM 求解一维饱和方程 - 吾爱随笔录

我想求解以下一维方程，以使用有限元模拟多孔介质中饱和度的时间演化（我是初学者）。渗透率和压力取决于饱和度。中的初始条件为： $K$ $p$ $s$ $\Omega = (0,l)$

s (t = 0) = 0.1, in Ω

$s(t=0)=0.1, \text{ in } \Omega$

s [0] = 0.5; left boundary

$s[0] = 0.5; \text{left boundary}$

\frac{\partial s}{\partial t} = - \frac{1}{ϕ} \nabla \cdot u u = \frac{K (s)}{μ} \nabla p (s)

$\frac{\partial s}{\partial t} = - \frac{1}{\phi} \nabla \cdot u \\ u = \frac{K(s)}{\mu} \nabla p(s)$

对于时间发展，我使用有限差分（显式欧拉）对时间导数进行离散化。我知道这不是最优的，因为需要求解线性系统才能获得新的时间步长，但我认为它应该适用于初步解决方案。

我使用了线性形状函数和标准的 Galerkin 方法（形状函数 = 权重函数）。我的近似解是：其中是线性形状函数。对于离散弱公式，指数表示时间，我得到：

s = \sum_{i} {\tilde{s}}_{i} N_{i}

$s = \sum_i \tilde s_i N_i$

N

$N$

n

$n$

M S^{n + 1} = M S^{n} + \frac{Δ t}{ϕ} \frac{K (S^{n})}{μ} D p^{n}

$MS^{n+1} = MS^{n} + \frac{\Delta t}{\phi}\frac{K(S^n)}{\mu}Dp^n$ 用质量矩阵：和刚度矩阵：。

m_{i j} = \int_{Ω} N_{i} N_{j} d Ω

$m_{ij} = \int_\Omega N_iN_jd\Omega$

d_{i j} = \int_{Ω} \nabla N_{i} \nabla N_{j} d Ω

$d_{ij} = \int_\Omega \nabla N_i \nabla N_j d\Omega$

然而，该方法似乎是错误的，因为解决方案没有收敛，并且在饱和跳跃处剧烈振荡。我查看了不同的问题，例如热方程，但发现它与我的系统不同。尤其是因为和的饱和依赖性是强非线性的（指数函数）。 $K$ $p$

会不会是标准 Galerkin 方法不适用？我将非常感谢任何提示。

PS：我可以使用基于粒子的模拟（SPH）成功模拟案例。