HbyA,phiHbyA,fvc::div(phiHbyA)计算错误问题

李东岳

Yes, 我觉得这个算例非常好，因为OKSS2上来就是高斯积分计算，也算是我对你这个感兴趣的动力，我会把这个放在课上分析。

由phiHbyA求和除以mesh.V()所得的结果，与输出的fvc::div(phiHbyA)相差了好几个量级！

用上面那个简单试试，可能你那个太复杂了，小数点很多，数也不工整。客观来讲，不会出现这个现象。

yhdthu

@李东岳前辈您说的这个算例我没异议，这个就和我的那个简单例子一样的，我目前已经掌握HbyA到fvc::div(phiHbyA)的计算方法。

但是，问题来了，我用相同的计算方法去计算我的第二个例子的某个cell，发现求和的时候出现了错误。

即，我手动计算的fvc::div(phiHbyA)为：-6.3169e-3
但是输出的结果却是1.90343e-6

如果只看phiHbyA求和，其结果是10^-10量级，cell的体积是10^-7量级，两个相除就是10^-3次方量级了；

输出的fvc::div(phiHbyA)的结果居然是10^-6量级，这可是除完体积之后的值，也就是说之前的phiHbyA的求和量级是10^-13量级

前辈，这就是我想表达的意思，为什么从phiHbyA就算不到fvc::div(phiHbyA)了？？

李东岳

所以就是说：你实际上已经知道手动计算的方式，但是你算的和代码结果不一样。如果可以的话，请把4号网格4个面的phiHbyA输出，我看看你是怎样 手动计算的fvc::div(phiHbyA) 的，就是说：把你计算-6.3169e-3的过程写出来。

输出的fvc::div(phiHbyA)的结果居然是10^-6量级

这个非常正常，匹配你的压力方程残差。
你可以看出来，这只不过是加减法问题，可能你手动计算存在一个非常小的bug。

yhdthu

@李东岳前辈，我的意思是：

我计算那个方形的4号cell全都没问题，以上第一个方形网格的算例给出了我所查找的所有参数；

有问题的是，用相同的计算方法，去计算第二个算例的47794号cell，从phiHbyA到fvc::div(phiHbyA)的计算出了问题。

我明白量级什么的无所谓，最重要的是数得算的对，但是现在手动算的和程序算的根本就对不上了，而且差了n个量级，这就说不过去了:crying:

yhdthu

@李东岳前辈具体的计算过程如下，这里只说有错误的第二个算例（cell_47794）

需要进行计算的场输出为：
HbyAVol = HbyA;
phiHbyASuf = phiHbyA;
sourcePoissionU = fvc::div(phiHbyA);

根据输出的phiHbyA，4个面的值如下：

E：0.000125002
W：0.000124669
N：-1.94136e-05
S：-1.908e-05
F 和 B 为empty，0

加和：E - W + N - S = 6e-10

mesh.V() = 9.49835e-8

fvc:div(phiHbyA) = (E - W + N - S) / mesh.V() = -0.0063169

但是，输出的fvc::div(phiHbyA) = 1.90343e-6

李东岳

E：0.000125002
W：0.000124669
N：-1.94136e-05
S：-1.908e-05

这个找的正确么

yhdthu

@李东岳前辈，这个百分百没问题，我是用paraview一个一个面点出来看的

李东岳

你是怎么在paraview查看face的值的？

yhdthu

@李东岳通过paraview找到cell编号，查找polymesh里面的owner和neighbor文件定位行号，打开surfacefield对应找就好啦～

李东岳

方法是对的，不过为什么错了？我也很奇怪。

yhdthu

@李东岳更奇怪的是，程序前半段计算phiHbyA都是对的，但是之后求div的时候（加法）却错了，难道加法也能算错:upset:

yhdthu

@李东岳前辈，解决了，原来是输出精度不够导致的，我的锅:lol:

CFD中文网