fvm::ddt(rho, T)与rho*fvm::ddt(T)的区别

[ChimneyFoam]

假设我又一个这样的方程：
$$\frac{\mathrm{\partial }(\rho T)}{\mathrm{\partial }t}=-\mathrm{\nabla }\cdot (\overrightarrow{U}\rho T)$$

其中$\rho$表示流体（比如水）的密度，根据状态方程及p,T计算出来。请问在离散的时候下面这两种写法的结果有什么区别？

1. fvm::ddt(rho, T) == -fvm::div(phi,T);

2. rho*fvm::ddt(T) == -fvm::div(phi,T);

[Samuel-Tu]

区别应该在于密度是否可压吧

[李东岳]

fvm::ddt(rho,T) = $\frac{\mathrm{\partial }\rho T}{\mathrm{\partial }t}$, rho*fvm::ddt(T) = $\rho \frac{\mathrm{\partial }T}{\mathrm{\partial }t}$

[李东岳]

我最近在更新rhoFoam系列求解器，做了更深的工作，发现ddt这一项还有点意思，假设$t$为当前时间步，欧拉格式：

• 在runTime++之前，fvc:ddt(T) = $\frac{T^t-T^{t-\mathrm{\Delta }t}}{\mathrm{\Delta }t}$;

• 在runTime++之后，如果T进行TEqn.solve()更新，fvc:ddt(T) = $\frac{T^{t+\mathrm{\Delta }t}-T^t}{\mathrm{\Delta }t}$

• 在runTime++之后，如果T没有进行TEqn.solve()之类的更新，fvc:ddt(T) = 0，因为$T^{t+\mathrm{\Delta }t}=T^t$, 同时T.oldTime()也是$T^t$。，

主要是代码中的vf，vf.oldTime()需要跟时间步对一下。

所以fvm::ddt(rho, T)还需要考虑rhoEqn是否更新：

• 如果在runTime++之后，进行了rhoEqn.solve()，那么fvm::ddt(rho, T) = $\frac{{\rho }^{t+\mathrm{\Delta }t}{T}^{t+\mathrm{\Delta }t}-{\rho }^t{T}^t}{\mathrm{\Delta }t}$

• 如果在runTime++之后，没有进行rhoEqn.solve()，那么fvm::ddt(rho, T) = $\frac{{\rho }^t{T}^{t+\mathrm{\Delta }t}-{\rho }^{t-\mathrm{\Delta }t}{T}^t}{\mathrm{\Delta }t}$

• 前者看起来更耦合在一起

二者都不同于rho*fvm::ddt(T)，因为其等于${\rho }^{t+\mathrm{\Delta }t}\frac{T^{t+\mathrm{\Delta }t}-T^t}{\mathrm{\Delta }t}$（进行了rhoEqn.solve()），或${\rho }^t\frac{T^{t+\mathrm{\Delta }t}-T^t}{\mathrm{\Delta }t}$（没进行rhoEqn.solve()）

[ChimneyFoam]

@东岳 太感谢了！您这次来了一个非常透彻的分析，我的疑问已全部被解决。

[闻久STU]

@李东岳 哇！！！透彻！！mark了，非常感谢

[Krasus]

感谢解读，学习到了@李东岳

[李东岳]

我在无痛苦ns方程笔记里面增加了对这一行代码的解释：

fvScalarMatrix alpha1Eqn
(
    fvm::ddt(alpha1) - fvc::ddt(alpha1) // WHY???
  - fvm::laplacian(alpha1alpha2f()*pPrimeByA_(), alpha1)
);

alpha1Eqn.relax();
alpha1Eqn.solve();

尤其是其中的

 fvm::ddt(alpha1) - fvc::ddt(alpha1) // WHY???