@李东岳 在 模拟风沙两相流，建模时带孔的网子怎样处理 中说：
@贾光普 后来这个工作有什么进展么
讲道理我现在在做类似于这个玩意的东西，只不过我做的是气液+颗粒三相的东西，现在卡在筛选颗粒这块好久了。目前已经确定筛网用多孔介质做了，nasa他们有很久很久之前的阻力系数测试的数据，但是这个颗粒筛选我在想要不要用liggghts重新写一个带筛选功能的反射边界。
1
16c
@16c
帖子
-
模拟风沙两相流，建模时带孔的网子怎样处理 -
codedsource源项不收敛然后又去pimpleFoam跑了一下，代码改成这样的：
momentumSource { type vectorCodedSource; active yes; name sourceTime; vectorCodedSourceCoeffs { selectionMode all; // cellZone pZone; fields (U); codeInclude #{ #}; codeCorrect #{ // Pout<< "**codeCorrect**" << endl; #}; codeAddSup #{ // Pout<< "**codeAddSup**" << endl; // const vectorField& C = mesh_.C(); const scalarField& V = mesh_.V(); vectorField& Usource = eqn.source(); const vectorField& U = mesh().lookupObject<volVectorField>("U"); // const scalarField& Rho = mesh().lookupObject<volScalarField>("rho"); const scalarField& nu = mesh().lookupObject<volScalarField>("nu"); // const scalarField& magU = mag(U); scalar A = 1e2; scalar B = 1e2; // vector C(0,1e4,0); forAll(V,i) { const scalar x = mesh_.C()[i][0]; const scalar y = mesh_.C()[i][1]; if(x < 0.5 && x > 0 && y < 0.5 && y > 0.45) { Usource[i] += (1e-5 * A + mag(U[i])* B * 0.5 ) * U[i]* V[i]; // Usource = (A * U[i] + B * mag(U[i]) * U[i]) * V[i]; // Usource[i] += - C * V[i]; } } // Info << "***codeAddSup***" << nl; #}; codeSetValue #{ // Pout<< "**codeSetValue**" << endl; #}; // Dummy entry. Make dependent on above to trigger recompilation code #{ $codeInclude $codeCorrect $codeAddSup $codeSetValue #}; } sourceTimeCoeffs { $vectorCodedSourceCoeffs; } }跑完的速度场是这样：
就完全处于A和B只能在1e2这个数量级，但凡再大一点就直接浮点溢出了。。。完全不知道是咋回事。。
-
codedsource源项不收敛interFoam跑完了的图是这样的
-
codedsource源项不收敛给源项公式改了一下，但是发现一个非常难受的事情，本身这个东西是仿照darcy-forchheimer写的，但是openfoam原生的模型的D和F可以提到1e10这么高的数量级，但是我这个写的感觉只能到1e2这样子，1e3开始就会浮点溢出，但是1e2这个水平只是能看出来这个特定区域确实有阻力存在，并不能实现极大多孔介质阻力并形成类固体区域的状态，所以就很好奇这个是什么情况，新改的代码是这样：
momentumSource { type vectorCodedSource; active yes; name sourceTime; vectorCodedSourceCoeffs { selectionMode all; // cellZone pZone; fields (U); codeInclude #{ #}; codeCorrect #{ // Pout<< "**codeCorrect**" << endl; #}; codeAddSup #{ // Pout<< "**codeAddSup**" << endl; // const vectorField& C = mesh_.C(); const scalarField& V = mesh_.V(); vectorField& Usource = eqn.source(); const vectorField& U = mesh().lookupObject<volVectorField>("U"); const scalarField& Rho = mesh().lookupObject<volScalarField>("rho"); const scalarField& nu = mesh().lookupObject<volScalarField>("nu"); // const scalarField& magU = mag(U); scalar A = 1e3; scalar B = 1e3; // vector C(0,1e4,0); forAll(V,i) { const scalar x = mesh_.C()[i][0]; const scalar y = mesh_.C()[i][1]; if(x < 0.5 && x > 0 && y < 0.5 && y > 0.45) { Usource[i] += (nu[i] * A + mag(U[i])* B * 0.5 ) * Rho[i] * U[i]* V[i]; // Usource = (A * U[i] + B * mag(U[i]) * U[i]) * V[i]; // Usource[i] += - C * V[i]; } } // Info << "***codeAddSup***" << nl; #}; codeSetValue #{ // Pout<< "**codeSetValue**" << endl; #}; // Dummy entry. Make dependent on above to trigger recompilation code #{ $codeInclude $codeCorrect $codeAddSup $codeSetValue #}; } sourceTimeCoeffs { $vectorCodedSourceCoeffs; } } -
openfoam在多孔介质是怎么实现仅有阻力的？最近在用codedsource写个阻力源项，但是写完之后感觉效果都不太好，昨天用paraview看了一眼矢量图，发现用多孔介质模型的区域是没有速度矢量的，但是用codedsource的话从t=0开始就存在速度矢量，然后导致介质满天飞。所以openfoam是怎么实现让多孔介质模型明明作为一个动量源项但是不沾任何速度矢量的。。。
-
codedsource源项不收敛@李东岳
这个东西我也很奇怪，就是我看fluid mechanics 101的多孔介质视频 给出的多孔介质公式是
这个形式的，然后看动量方程里面各项基本都是力的形式，我之前就想当然在后面的公式里面多乘了个体积V,后来看of帮助文档里面写的是
我就又改成上面那个形式，但是感觉好像写完了之后无论哪个形式也不收敛。。。 -
codedsource源项不收敛我在fvoptions里面想设置一个区域性的阻力源项 ，但是因为需要实现的结构比较特殊，所以没法用toposet去进行捕捉获取，然后知乎看见一个大佬用codedsource写的一个特定区域的动量源项合计自己试试去写一个这种，为了尝试设置的阻力源项可不可用我就toposet了一个区域合计试试，但是无论我的源项是纯数字还是仿darcy-forchheimer公式类型的他用interFoam跑都不收敛，找了好几天论坛视频啥的也找不到原因，vscode报错也只是给代码报错，偶尔代码看起来没错误但是跑起来就是不收敛。现在目前写的是如下这样的，因为回头要耦合dem，所以of版本是5.x。
momentumSource { type vectorCodedSource; active yes; name sourceTime; vectorCodedSourceCoeffs { selectionMode cellZone; // all; cellZone pZone; fields (U); codeInclude #{ #}; codeCorrect #{ // Pout<< "**codeCorrect**" << endl; #}; codeAddSup #{ // Pout<< "**codeAddSup**" << endl; // const vectorField& C = mesh_.C(); const scalarField& V = mesh_.V(); vectorField& Usource = eqn.source(); const vectorField& U = mesh().lookupObject<volVectorField>("U"); const scalarField& Rho = mesh().lookupObject<volScalarField>("rho"); const scalarField& nu = mesh().lookupObject<volScalarField>("nu"); const scalarField& magU = mag(U); scalar A = 10; scalar B = 10; forAll(V,i) { Usource -= (A * nu * Rho + B * Rho * magU / 2) * U ; } // Info << "***codeAddSup***" << nl; #}; codeSetValue #{ // Pout<< "**codeSetValue**" << endl; #}; // Dummy entry. Make dependent on above to trigger recompilation code #{ $codeInclude $codeCorrect $codeAddSup $codeSetValue #}; } sourceTimeCoeffs { $vectorCodedSourceCoeffs; } }请各位大佬各位巨佬帮着指点一下哪里出现了问题。
-
interfoam计算出现了奇怪的结果，猜测是压力的问题，但是不知道具体改怎么做。。。@李东岳 我因为可能要做cfdem所以用的5.x跑的
5.x改了好多次网格尺寸和dt都在浮点溢出，后来看网上教程给piso改成pimple模式 倒是不溢出了，但是这个结果就很离谱了
-
interfoam计算出现了奇怪的结果，猜测是压力的问题，但是不知道具体改怎么做。。。 -
interfoam计算出现了奇怪的结果，猜测是压力的问题，但是不知道具体改怎么做。。。在openfoam里用interFoam算了一个很简单的cas，大概就是一个平板前面有个inlet,往里面注水，想看这个水流的情况，但是出现了非常离谱的状况，我先放上结果图
我设置的速度方向是(1 0 0),因此我觉得他不应该出现向上喷射的这个状态。
同样的我又在fluent里面简单跑了一下这个cas，边界条件设置都是相同的，结果是这样的：
今天着重的了解了一下这个p_rgh文件，大概可以理解为引入的一个方便求解的中间参量，等于原来的p-rhogh。
两者设置里面最唯一出现的差别可能就是压力，fluent的设置更加傻瓜式一点，直接对边界区域的相对压力给定值就好了，这个p_rgh原则上感觉不需要去关心具体是多少，但是总感觉这个东西好像跟不同的坐标有一定的关系，比方说坐标（ 1 1 0）的点和坐标（ 1 1 1 ）的点在初场设置时想让他们的静压或者说P文件中的压力为0的时候，需要设置在的p_rgh中的值应该是不一样的。
感觉出现这个问题也可能是我openfoam这边设置的问题。所以也得请求各位大佬帮着看一下，openfoam的设置文件我放在下面了。
test2.zip因为文件稍微有点大我就给网格删了，不过都是均匀网格所以画起来也不麻烦，直接blockMesh就好，不过需要topoSet+createPatch一下来创造一个inlet来着
-
有关p_rgh我自己理解这个这个p_rgh应该是其他求解器0文件下的p-rg·h，但是interFoam似乎只能读取p_rgh而不能读取p文件，但是在对边界和整个初场只已知p的时候p_rgh怎么设置，尤其是inlet和outlet这种在某些时候需要给定一些数值的时候。
-
openfoam对网格密度的敏感度喝适配度是什么样的最近在跑几个非常简单的case，发现网格画的比较精细的时候反而经常库朗数爆炸，网格画的相对粗糙的时候反而跑的杠杠的。 是与电脑配置有关还是说求解器有个适配的敏感度啥的，我是一点头绪也没有。
而且报错之后浮点溢出的东西在我看起来就像是一坨乱码，unset FOAM_SIGFPE之后也没有具体的错误原因，这个报错怎么看有没有啥关键字或者资料啥的，还得请各位大佬提点提点。 -
toposet之后出现空集。。。@学流体的小明 在 toposet之后出现空集。。。 中说：
感觉是topoSet选择网格，提取各种面的时候有些问题了，我没用过这方面的功能，不太了解。
你如果只是想把边界面提取出来的话，在blockMesh当中就可以设置吧。@李东岳 在 toposet之后出现空集。。。 中说：
豁牙子注意一下toposet点位置的精度。
感谢两位大佬提点，豁牙子问题解决了...
-
toposet之后出现空集。。。@学流体的小明 在 toposet之后出现空集。。。 中说：
你blockMesh文件convertToMeters 0.001;
实际就是把下面的坐标乘以0.001画网格，画出来的计算域是(0,0,0)到(2,0.5,0.05)。
而你toposet直接从(0,0,10)开始画box，当然找不到了。大哥 按你说的改了topoSetDric , box的坐标和网格尺寸我也确定了是能对应上的，但是createPatch并ParaFoam之后出现了inlet这个面变成了豁牙子。。。
改之后的topoSet是这样的：
/*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 5 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; location "system"; object topoSetDict; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // actions ( { name porous; type cellSet; action new; source boxToCell; sourceInfo { box (0 0 10e-3) (2000e-3 500e-3 12e-3); } } { name porousZone; type cellZoneSet; action new; source setToCellZone; sourceInfo { set porous; } } { name inletA; type cellSet; action new; source boxToCell; sourceInfo { box (0 0 12e-3) (4e-3 500e-3 22e-3); } } { name inletB; type faceSet; action new; source cellToFace; sourceInfo { set inletA; option all; } } { name inletB; type faceSet; action subset; source patchToFace; sourceInfo { name walls; } } { name outletA; type cellSet; action new; source boxToCell; sourceInfo { box (1.996 0 12e-3) (2 500e-3 50e-3); } } { name outletB; type faceSet; action new; source cellToFace; sourceInfo { set outletA; option all; } } { name outletB; type faceSet; action subset; source patchToFace; sourceInfo { name walls; } } ); // ************************************************************************* //createPatch是这样的：
/*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 5 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; object createPatchDict; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // Do a synchronisation of coupled points after creation of any patches. // Note: this does not work with points that are on multiple coupled patches // with transformations (i.e. cyclics). pointSync false; // writeCyclicMatch false; // Patches to create. patches ( { name inlet; patchInfo { type patch; } constructFrom set; set inletB; } { name outlet; patchInfo { type patch; } constructFrom set; set outletB; } ); // ************************************************************************* //这个豁牙子问题我之前跟着阿B上田东的视频做案例的时候也出现了 不知道是什么玩意 而且case跑完之后他明显是个实际存在的数值上的豁牙子，不是paraview显示问题
-
toposet之后出现空集。。。toposet之后polymesh下set文件夹的文件内全部显示是空集，终端显示识别了你的范围是XXX但是这个范围内狗der没有，我是真不到哪里出问题了。
toposetdric代码是这样的：/*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 5 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; location "system"; object topoSetDict; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // actions ( { name porous; type cellSet; action new; source boxToCell; sourceInfo { box (0 0 10) (2000 500 15); } } /* { name porousZone; type cellZoneSet; action new; source setToCellZone; sourceInfo { set porous; } } */ { name inletA; type cellSet; action new; source boxToCell; sourceInfo { box (0 0 15) (0.1 500 25); } } { name inletB; type faceSet; action new; source cellToFace; sourceInfo { set inletA; option all; } } { name inletB; type faceSet; action subset; source patchToFace; sourceInfo { name walls; } } { name outletA; type cellSet; action new; source boxToCell; sourceInfo { box (1999.9 0 15) (2000 500 50); } } { name outletB; type faceSet; action new; source cellToFace; sourceInfo { set outletA; option all; } } { name outletB; type faceSet; action subset; source patchToFace; sourceInfo { name walls; } } ); // ************************************************************************* //终端是这样说的：
openfoam@openfoam-virtual-machine:~/OpenFOAM/openfoam-5.x/run/1$ topoSet /*---------------------------------------------------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 5.x | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ Build : 5.x-7f7d351b741b Exec : topoSet Date : Nov 25 2024 Time : 09:45:57 Host : "openfoam-virtual-machine" PID : 181025 I/O : uncollated Case : /home/openfoam/OpenFOAM/openfoam-5.x/run/1 nProcs : 1 sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE). fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 10) allowSystemOperations : Allowing user-supplied system call operations // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // Create time Create polyMesh for time = 0 Reading topoSetDict Time = 0 mesh not changed. Created cellZoneSet porous Applying source boxToCell Adding cells with center within boxes 1((0 0 10) (2000 500 15)) cellZoneSet porous now size 0 Created cellSet inletA Applying source boxToCell Adding cells with center within boxes 1((0 0 15) (0.1 500 25)) cellSet inletA now size 0 Created faceSet inletB Applying source cellToFace Adding faces according to cellSet inletA ... faceSet inletB now size 0 Read set faceSet inletB with size 0 Applying source patchToFace Adding all faces of patch walls ... Found matching patch walls with 31250 faces. faceSet inletB now size 0 Created cellSet outletA Applying source boxToCell Adding cells with center within boxes 1((1999.9 0 15) (2000 500 50)) cellSet outletA now size 0 Created faceSet outletB Applying source cellToFace Adding faces according to cellSet outletA ... faceSet outletB now size 0 Read set faceSet outletB with size 0 Applying source patchToFace Adding all faces of patch walls ... Found matching patch walls with 31250 faces. faceSet outletB now size 0 EndblockMesh是这样的：
/*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 5 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | FoamFile { version 2.0; format ascii; class dictionary; object blockMeshDict; } convertToMeters 0.001; vertices ( (0 0 0) // 0 (2000 0 0) // 1 (2000 500 0) // 2 (0 500 0) // 3 /* (0 0 15) // 4 (2000 0 15) // 5 (2000 500 15) // 6 (0 500 15) // 7 */ (0 0 50) // 4 (2000 0 50) // 5 (2000 500 50) // 6 (0 500 50) // 7 // (0 0 25) // 12 // (0 500 25) // 13 ); blocks ( hex (0 1 2 3 4 5 6 7) (500 125 25) simpleGrading (1 1 1) // hex (4 5 6 7 8 9 10 11) (500 125 35) simpleGrading (1 1 1) ); edges ( ); boundary ( walls { type wall; faces ( (0 1 5 4) (2 3 7 6) (0 4 7 3) (1 2 6 5) ); } atmosphere { type patch; faces ( (4 5 6 7) (0 1 2 3) ); } mergePatchPairs ( );我搞这个破玩意好几天了，一点头绪没有，还得求着各位大佬帮帮忙瞅一眼。。。。