CFD中文网

@Micro final就是PIMPLE循环的最后一步。

@hongjiewang 看起来应该没问题，中间的网格密是因为用refinement region加密过了。

@hongjiewang 不需要对两个区域分别划分，利用 sHM 可以直接同时生成一整套网格。

@hongjiewang 不需要merge，只需要一个 case 文件夹就行。

@hongjiewang 在 不规则多孔介质的网格划分 中说：

对于一些简单几何模型直接采用topoSet是可行的，对于不规则的还是需要借助stl进行区域的设定。

你在另一个帖子里问我了，我在这里回复你吧。

在 snappyHexMeshDict 里面可以通过 STL 文件（必须封闭）指定 cellZone。例如：

geometry
{
    sphere
    {
        type triSurfaceMesh;
        file "sphere.stl";
    }
};

refinementSurfaces
{
    sphere
    {
        level (0 0);           
        faceZone sphere1;
        cellZone sphere1;
        cellZoneInside inside;
    }
}

@Stan同学 在 关于snappyHexMesh 中说：

@wwzhao
您好，想向您请教个小问题。

前面这个网址 https://sites.google.com/site/snappywiki/ 里讲到如果要添加边界层，需要指定具体的patch，但是我不太清楚stl文件内部是如何定义各个patch的，请问有什么办法或者用什么软件来打开stl文件从而获取里面的patch信息吗？

在生成stl文件时应该给solid命名，以区分stl中的不同patch。

@Stan同学 在 关于snappyHexMesh 中说：

另外想请问您说的每个patch单独保存一个stl文件，这个用CAD软件可以实现吗 ？
（之前一直是用的CAD导出stl文件，应该是您说的整个区域为一个stl文件）

这个可以通过FreeCAD实现。

插值只对patch进行，和是否并行无关。并行对local patch进行插值。

@Samuel-Tu primitivePatchInterpolation

faceProcAddressing是并行用的，将local face id转换成global face id。

@heike256 试试 spack。

写个类调用下fftw的函数接口就好了

@gyzhangqm

很精彩的推导！解开了我一直以来的疑惑。

不过有一点我不认同：

这里有个trick，$a _ { N } = \dot { m } _ { f } ( 1 - \varpi ) _ { f }$，为什么对于编号大于Owner的邻居单元，upper 存的是$\dot { m } _ { f } ( 1 - \varpi ) _ { f } $，而对有编号小于Owner的邻居单元，lower存的是 $- \varpi _ { f } \dot { m } _ { f } $呢?这是因为计算面心的系数$\varpi$时在总是认为面左侧是Owner，右侧是Neighbour.

这里的trick（正负号）应该是求flux的时候设置的，而不是求weight。

我换个说法理一下思路：

对于单元 $P$，有 $k$ 个面，也就有 $k$ 个邻居单元，记为 $N_1$ ... $N_k$。

每个面心的量通过面的 owner 和 neighbour 单元体心插值得到（注意这里的 owner、neighour 和上面的 $P$、$N$ 完全不相关，为了避免歧义，上面的单元用 $P$ 而不是用 $O$ 表示）

$$
\phi_ f = \varpi \phi_{own} + (1 - \varpi) \phi_{nei}
$$

于是对流项的积分变为

$$
\sum_{f=1}^{k} (\rho \mathbf { u } \phi) _ {f} \cdot \mathbf S_f = \sum_{f=1}^k \dot m_f \phi_f
$$

上式中的 $\dot m$ 为质量通量，这个是由速度点乘面法向得到的。具体操作为 fvc::flux(U)，也就是 mesh.Sf() 和 U 的点乘。Sf() 是一个矢量场，它的方向很有讲究。它被规定为从 own 指向 nei。而高斯定理的 $\mathbf S_f$ 规定方向指向封闭单元外。这两个定义是由矛盾的。因此需要上面的trick来调整正负号。

先看非对角元素。遍历 P 的所有面单元：

• 若该面单元的 own 为 P， nei 为 N，Sf() 与 $\mathbf S_f$ 同方向
  这种情况下，N 的系数在 upper，因此 upper 需要加上 $\dot m(1 - \varpi)$。同时 P 的主对角元素加上 $\dot m \varpi$。
• 若该面单元的 nei 为 P，own 为 N，Sf() 与 $\mathbf S_f$ 反方向
  这种情况下，N 的系数在 lower，原本 lower 需要加上 $\dot m \varpi$，调整负号后，最终加上的是 $-\dot m \varpi$。同时 P 的主对角元素加上 $-\dot m(1 - \varpi)$。

以上过程对应以下两行代码：

    fvm.lower() = -weights.primitiveField()*faceFlux.primitiveField();
    fvm.upper() = fvm.lower() + faceFlux.primitiveField();

再看主对角元素。

对于邻居单元 N ，在计算其对应行的非对角系数时，所用插值公式和 P 相同，即还是

$$
\phi_ f = \varpi \phi_{own} + (1 - \varpi) \phi_{nei}
$$

也就是说以上公式要用两次，分别是在 P 和 N 作为主对角元素的行上使用。

在循环到 N 与 P 共享的面单元时：

• 若该面单元的 own 为 P， nei 为 N，Sf()方向为 P->N ，$\mathbf S_f$ 方向为 N->P，二者反方向
  这种情况下，P 的系数为 $-\dot m \varpi$。
• 若该面单元的 nei 为 P，own 为 N，Sf() 方向为 N->P， $\mathbf S_f$ 方向也为 N->P，二者同方向
  这种情况下，P 的系数为 $\dot m(1 - \varpi)$ 。

注意以上两种情况，P 作为 N 的非对角系数和作为其本身的主对角系数正好是相反数，且处于同一列上。由于遍历一次面单元将所有网格单元（包括P和其所有邻居单元）的非对角元素都计算好了，为了避免重复计算，所以直接用该列的非对角元素加和并取相反数得到主对角元素的值。

即对应代码：

    fvm.negSumDiag();

也就是说 negSumDiag 计算的确实是列，而不是行。我之前会由系数矩阵对称的错觉是因为采用的验证算例都是采用线性插值计算面心物理量，这使得插值系数 $\varpi=0.5$，正好造成系数矩阵反对称的假象。

@gyzhangqm

$$
a_C=-\sum_{f:nb(C)}{a_N}+\sum_{f:nb(C)}{\dot m_f}
$$

上面这个式子是不是漏了 $\phi_C$ 和 $\phi_N$？$\dot m_f$ 是不是指的 Ax=b 中的 b？

我的理解应该是

$$
a_C\phi_C = - \sum_{f:nb(C)}a_N\phi_N + b_C
$$

以上公式为计入边界条件后的表达式， $b_C$ 不为 0，且 $a_C$ 的值被边界条件影响。

negSumDiag 只对 diag， lower 和 upper 做了操作，并没有计入边界条件的影响，所以其过程并不对应以上公式。

y+>5时wallGradU计算出来的并不是真实的法向速度梯度，这种情况下通过wallGradU得到的壁面切应力也是不正确的。

@Samuel-Tu 在 编译失败求助 中说：

void crtMyNu(const labelUList& myList,scalar myNu)

这句改成 void singlePhaseTransportModel::crtMyNu(const labelUList& myList,scalar myNu) 试试？

probe 需要获取物体的运动信息才能确定如何移动，所以光改这个还不够。

@七辆战车

如何在dynamicMeshDict中实现用弹簧连接双浮体呢，找了找网上的资料，没能够找着相关资料

这个应该没有现成模块，需要自己开发。

只有floatingObject2可以动，floatingObject1不能动

sixDoFRigidBodyMotion 不支持多物体运动，你可以试试 rigidBodyMotion。

1. 找不到centreOfMass肯定是字典文件错了，检查字典文件

2. linearSpring 没法实现双浮体相连的功能

@VamosZhi simpleFoam 系列求解器都是稳态计算，oversetSimpleFoam 是不是不支持动网格？

@VamosZhi 这个和overset没有关系，应该是别的原因导致的