额，我应该早点看到这个帖子的，今天面试就问了这个问题，我有点懵，答得不太好。。。

@李东岳 感谢李老师！

原来是这样，感觉很神奇啊。请教李老师，为什么要对源项分开进行显式隐式的离散？

对对对，这样是对的，我之前把应变率和湍动能搞混了，谢谢李老师了！

李老师的“CFD中的大涡模拟”中的公式19，\begin{equation}
\nu_{\mathrm{SGS}}=\rho C_{\mathrm{SGS}} \Delta \sqrt{k_\mathrm{SGS}}
\label{Smag}
\end{equation}
CSGS和delta上面是不是少了平方？看了两本书都是有平方的。还有运动粘度的话是不是多乘了一个密度?

@Cr 原来是在这里啊，感谢感谢！

在各种源项类的声明和实现中经常能看到addSup函数，功能为将源项加入到控制方程中，如在actuationDiskSource.H中：

 //- Source term to momentum equation
virtual void addSup
 (
         fvMatrix<vector>& eqn,
         const label fieldi
 );

那么请问这个addSup函数是在什么时候被调用的？找了很久没有找到。在网上查了一些资料，基本可以确定这个函数是在求解器形成U方程的时候被调用了，如在simpleFoam的UEqn.H中：

     tmp<fvVectorMatrix> tUEqn
     (
         fvm::div(phi, U)
       + MRF.DDt(U)
       + turbulence->divDevSigma(U)
      ==
         fvOptions(U)
     );

大家都说fvOptions(U)就用到了addSup函数，可是这里addSup函数是如何被调用的？希望有大神能帮忙解答！

时隔九个月，今天重新看这个帖子，终于搞清楚了，记录一下。
确实是像李老师说的那样，Usource()就是源项，这里是手动加的源项，之前一直以为是系数矩阵的非对角元。
Udiag()代表U方程的对角线元素，来自IduMatrix，而Usource()是源项矩阵，来自fvMatrix。上面的代码之所以会出现Udiag()，其实就是像李老师说的那样，是为了增加对角占优。而多孔隙模型的源项就是从Usource加入的：

Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);

这里Cd就是公式（2）中括号中的部分。

@秋山结月
1.对的，最好另外再建立新的几何体，用这个新的几何体包裹要refine的局部区域。
2.这个问题倒是没有遇到过，是不是stl文件的问题？
3.snap结束后不应该直接没掉啊，即便是没有需要施加边界层的区域，也会有“No layers to generate ...”，这样的信息。看看你后面有没有这样的提示信息？

个人的一个经验是善用refinementRegions关键字，在需要加密的地方建立refinementRegions，将level提高些，应该可以优化一下网格。还有一个地方是可以试试调整surfaceFeatureDict中的includedAngle，可以优化角度变化比较剧烈处的网格，比如你第一张图那里。另外可以适当将maxLocalCells和maxGlobalCells提高些。这个还是熟能生巧，多试试改变那些参数。
SHM是分步来处理网格的，先是snap，最后才是加边界层。所以边界层网格在最后一个时间步内，不知道你有没有查看最后一个时间步的网格。