关于cyclic boundary condition
-
能不能帮我看一下这个问题?谢谢~
http://www.cfd-china.com/topic/1150/openfoam中cyclic周期性边界的问题 -
关于cyclic boundary最近有一些发现,跟楼主和各位高手探讨一下,看看能否有所启发。
本人的测试用例为大佬李东岳翻译的openfoam用户指南-9里面的顶盖驱动流,差别是左右边界设置为cyclic边界,设置方法为陈与论帖子(https://zhuanlan.zhihu.com/p/84137342 )前半段那个简单的方法,网格x和y两个方向各画了3个网格。openfoam版本用的是OpenFOAM-9.
针对这个案例,在研究对流项的离散时,在gaussConvectionScheme.C中fvmDiv函数里发现,计算边界对矩阵A和源项的影响时,fvm.internalCoeffs()[patchi]和fvm.boundaryCoeffs()[patchi]的计算没有区分边界是否为coupled(),是否coupled()是通过调用不同类型的valueInternalCoeffs(pw)和valueBoundaryCoeffs(pw)实现的。通过valueInternalCoeffs(pw)和valueBoundaryCoeffs(pw)的返回值我们可以看到,对于cyclic的两个边界,这俩函数在coupledFvPatchField.C中计算,分别返回了
Type(pTraits<Type>::one)*w; Type(pTraits<Type>::one)*(1.0 - w);
由于
fvm.internalCoeffs()[patchi] = patchFlux*psf.valueInternalCoeffs(pw); fvm.boundaryCoeffs()[patchi] = -patchFlux*psf.valueBoundaryCoeffs(pw);
结合phi_f = A*phi_c+B(https://zhuanlan.zhihu.com/p/609043262 )的理论,相当于在cyclic边界上fvm计算面心的值时,只用到了边界单元的体心值phi_c,而没有用到周期边界另一头的单元值。分析其原因,这可能是由于openfoam矩阵存储方式的限制不得已而为之,因为非对角元是由lduAddressing来定位的,对于本案例而言,非对角上的非零元个数与内部面数目相同,为12个,如果在cyclic边界上计算面心的值时用到两个体心的值,则会新增一个非零元(相当于这个面变成了内部面),那么就出现了13个内部面,没办法在当前体系下存储。openfoam在处理周期边界fvm对流项的离散时相当于有一个窟窿。
那么,既然没有用到cyclic另一头的单元,OpenFOAM怎么实现周期边界的呢。在surfaceInterpolationScheme.C:287可以发现,在fvc类型计算插值时,对于coupled()的情况进行了分别处理,用到了patchNeighbourField(),相当于在计算通量时补上了前面的窟窿。但是感觉影响应该还是会有,因此我建议楼主考察一下,遇到的问题是否与openfoam这个窟窿有关。
最后,想跟大家探讨一下,假如想补上这个窟窿,有没有什么办法能够存储由于周期边界带来的矩阵中多的这个非零元素?这个问题非常关键,一是可以解决周期边界现有的问题,二是有可能解除当前openfoam面心插值时只有owner和neighbour两个体心值可用的限制。以上是个人的一些粗浅理解,希望大佬们指正!