CFD中文网

@派大星huangxiaojie

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  9
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/

forces1-gas
{
    type          forces;

    libs          ("libforces.so");

    writeControl  timeStep;
    timeInterval  1;

    log           yes;
    patches     (impeller_66 impeller_67);
    rho           rhoInf;     // Indicates incompressible
    log           true;
    rhoInf        700;          // Redundant for incompressible
    phase    gas;
     p p;
    CofR          (0 0.35 0);    // Rotation around centre line of propeller
    pitchAxis     (0 5 0);
}

// ************************************************************************* //

把每个相的扭矩值都分别处理，然后相加即可，跟实验值较为吻合

@李东岳 好的，谢谢老师

@李东岳 另外老师还想请教个问题，scalarTransportFoam中库朗数还需要小于1吗，还是说时间步长可以大一些让库朗数大一些

@李东岳 我的意思是把multiphaseEulerFoam算出来的U.particles和phi.particles一起复制过来作为scalarTransportFoam的初值，然后把原先phi.particles中的边界条件修改一下就可以了

@李东岳 老师，问题解决了，要把固体的通量场也作为初值给scalarTransportFoam，并且要把phi里面得fixedValue类型改为calculated就可以了，谢谢老师

@李东岳 或者说是不是多相流中的单相速度场并不能直接拿到scalarTransportFoam中使用，因为单相速度场没有相含率信息，例如一个网格的固含率为0，但是其还是会有一个固相速度，然而scalarTransportFoam仅读取速度场无法知道这一点，会不会导致质量不守恒的问题。

@李东岳 老师，我继续计算了4个小时也还是没有变化，出口一直没有浓度，而且T的迭代次数已经都是0了，就是说已经稳态了，可是这样如何得到停留时间分布？

@李东岳 老师，因为我想固体颗粒在液体中应该没有扩散，所以我只考虑了对流，扩散系数设了0，是这个问题吗

@李东岳 进口固定值1，其他wall和出口都是零梯度

并且算着算着迭代次数逐渐变为0了

首先使用multiphaseEulerFoam模拟了一个气液固三相连续搅拌釜，想获得固体相在釜中的停留时间分布，所以想将固体的速度场和通量场给到scalarTransportFoam中，从入口进入示踪剂以获得固体停留时间分布。这时遇到一个问题，整个釜的空时大概4个小时，但是计算大概1000s之后整个釜内的示踪剂浓度分布就不变化了，出口的浓度也达到稳定并且极低，请问这是什么问题？附上我的设置

ddtSchemes
{
    default         Euler;
}

gradSchemes
{
    default         Gauss linear;
}

divSchemes
{
    default         none;
    div(phi,T)     Gauss linearUpwind grad(T);
}

laplacianSchemes
{
    default         none;
    laplacian(DT,T) Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}

solvers
{
    T
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 1;
}
relaxationFactors
{

    equations
    {
        ".*"         1;
    }
}

@tidedrinker 你好，我发现of9中可以在force文件中单独定义相名称来提取相分数和对应的密度，输出的力和力矩如下：

forces forces write:
    sum of forces:
        pressure : (-10.813 0 -5.80003)
        viscous  : (-0.000457782 -0.000707157 0.000145593)
        porous   : (0 0 0)
    sum of moments:
        pressure : (-0.00441956 4.53963 -0.00703327)
        viscous  : (-3.21179e-05 7.86609e-05 -2.42533e-05)
        porous   : (0 0 0)

我的搅拌轴是y轴，所以我的扭矩是否就是4.53963N·M，但是这和实验值相差甚远，实验值大概只有0.06N·M，这是为什么？是不是还需要对这个力矩进行下一步处理才可以得到扭矩？需要除以密度嘛？

@tidedrinker 但是还有个问题，我在已经有了计算结果之后使用这个force进行后处理 -postProcess时，会报错说没有找到U和p，并不会自动读取U.*和p_rgh，请问您是怎么在多相流中使用的？

Reading field p_rgh

--> FOAM Warning :
    From function int main(int, char**)
    in file ./createFields.H at line 34
    Reading "/mnt/e/Desktop/limited-zeroG/constant/phaseProperties" from line 16 to line 144
    Pressure limits, pMin and pMax, are now read from "/mnt/e/Desktop/limited-zeroG/system/fvSolution/PIMPLE"
forces forces:
    Not including porosity effects


--> FOAM FATAL ERROR:
Could not find U, p

    From function void Foam::functionObjects::forces::initialise()
    in file forces/forces.C at line 205.

FOAM exiting

@tidedrinker 感谢，我现在在用of9的multiphaseEulerFoam，也就是说rho和rhoinf随便填一个数都没有影响是吗

各位大佬好，最近在做多相搅拌釜，想要提取桨叶的扭矩来计算搅拌功率。想要使用function中的force功能来提取桨叶上的扭矩，但是发现这个功能好像只适用于单相流。因为force文件需要输入rho，但是在多相流中由于相分布并不均匀且在瞬态变化，无法获得桨叶附近的rho值，所以该方法是否无法计算多相流的扭矩？

forces
{
    type          forces;

    libs          ("libforces.so");

    writeControl  timeStep;
    timeInterval  1;

    log           yes;

    patches       ("impeller*");
    rho           rhoInf;     // Indicates incompressible
    log           true;
    rhoInf        1;          // Redundant for incompressible

    CofR          (0 0 0);    // Rotation around centre line of propeller
    pitchAxis     (0 1 0);
}

另外拜读李东岳老师的博士论文发现李老师在论文中提取出了液-液两相流中的搅拌扭矩，说明多相流的扭矩是可以提取的，但是具体如何实现在论文中并未赘述，能否有大神帮忙解答，感谢！

@15735420621 可以参考自带案例mixerVessel2D，旋转域使用topoSet进行设置，旋转参数在MRFProperties中进行设置

@知行 ok，我的网格都是ICEM画的，我加密试试，谢谢

@知行 感谢耐心解答！首先网格的质量如下：

Mesh stats
    points:           182724
    faces:            532288
    internal faces:   516536
    cells:            174804
    faces per cell:   6
    boundary patches: 3
    point zones:      0
    face zones:       1
    cell zones:       3

Overall number of cells of each type:
    hexahedra:     174804
    prisms:        0
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    0
    polyhedra:     0

Checking topology...
    Boundary definition OK.
    Cell to face addressing OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
    Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces...
    Patch               Faces    Points   Surface topology
    walls               10840    10915    ok (non-closed singly connected)
    impeller1           1324     1326     ok (closed singly connected)
    outlet              3588     3663     ok (non-closed singly connected)

Checking geometry...
    Overall domain bounding box (-0.109995 -1.11019e-16 -0.109997) (0.109995 0.22 0.109997)
    Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
    Mesh has 3 solution (non-empty) directions (1 1 1)
    Boundary openness (5.83486e-16 -2.30396e-16 4.46812e-17) OK.
    Max cell openness = 3.29007e-16 OK.
    Max aspect ratio = 5.02531 OK.
    Minimum face area = 1e-06. Maximum face area = 2.73916e-05.  Face area magnitudes OK.
    Min volume = 4.64288e-09. Max volume = 1.25165e-07.  Total volume = 0.0083523.  Cell volumes OK.
    Mesh non-orthogonality Max: 39.7405 average: 7.14766
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.506949 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End

最大非正交性为39，平均为7，正交性还不错。
其次，几何是一个直径和高都是22cm的釜，桨叶直径8cm。
关于网格数量的话，由于目前还在案例调试阶段，所以先用17w的网格算起来比较快，五六十万的网格我也算过，也会出现这种现象。LES的计算代价有点大，所以还是想考虑RANS。
最后，在所有壁面处我使用了湍流模型kEpsilon的壁面函数，所以是不是对壁面处的网格厚度依赖性不是很大？我尝试把壁面加密一下试试。谢谢！

@李东岳 老师，MRF区域应该没啥问题，沿Y轴旋转，650rmp。
然后速度场看上去也是在转的：