[OlegSutyrin]

Hi, guys!
I aim to simulate inviscid compressible supersonic (M=5) flow over some complex solid axisymmetric body with rhoCentralFoam v2.3.1. I have some problems with it, so now I test some simpler cases instead: <1> cylindrical analogue of forwardStep (see "stepCyl.jpg") and <2> the same with equal backwardStep behind it ("step2Cyl.jpg"). Meshes are made with gmsh.

The grid is 1-layered wedge of 5 degrees wide, the bottom edge is the symmetry axis (not exactly, r=0.001 there). The air flows from left to right in the left part of the flow (U=745 m/s, see rough image of it on "initial.jpg") and is at rest everywhere else. Pressure (~850 Pa) and temperature (~75 K) is uniform. Inlet (left) is fixedValue patch; top and outlet (right) are zeroGradient patches; all others are walls with zeroGrad for P, T and slip for U.

The problem is that first case works perfectly well (shockwave generated by velocity gap passes to the right and flow eventually becomes stationary in 2-3 milliseconds), but the second suddenly one breaks down shortly after 0.05 milliseconds (everything becomes "nan"). The only difference (apart from extremely slight mesh differences) is the backward step added in the second case, and it is situated in the area that is completely at rest! At 0.05-0.06msec shockwave is still far to the left of it.

I'm completely lost, what may be the problem here? The only thing I know that pressure suddenly becomes negative at a single cell or face, but I don't know which one - can't see it in the paraView or Tecplot.

If you want to test it, cases and meshes may be downloaded here: tests.zip. "1.lay" files are layouts for Tecplot 360, in case you prefer it.

[OlegSutyrin]

Update:
I've redone the mesh with blockMeshDisct as follows:

PHP Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.1                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

convertToMeters 1;

vertices
(
    (-0.01 0 0) //0
    (0.1 0 0) //1
    (0.1 0.04642577085690893 -0.002026992930867164) //2
    (0.1 0.04642577085690893 0.002026992930867164) //3
    (-0.01 0.04642577085690893 -0.002026992930867164) //4
    (-0.01 0.04642577085690893 0.002026992930867164) //5
    (0.1 0.07992385772654863 -0.00348955098922688) //6
    (0.1 0.07992385772654863 0.00348955098922688) //7
    (-0.01 0.07992385772654863 -0.00348955098922688) //8
    (-0.01 0.07992385772654863 0.00348955098922688) //9
    (0.25 0.04642577085690893 -0.002026992930867164) //10
    (0.25 0.04642577085690893 0.002026992930867164) //11
    (0.25 0.07992385772654863 -0.00348955098922688) //12
    (0.25 0.07992385772654863 0.00348955098922688) //13
    (0.25 0 0) //14
    (0.4 0 0) //15
    (0.4 0.04642577085690893 -0.002026992930867164) //16
    (0.4 0.04642577085690893 0.002026992930867164) //17
    (0.4 0.07992385772654863 -0.00348955098922688) //18
    (0.4 0.07992385772654863 0.00348955098922688) //19
);

blocks
(
    hex (0 1 2 4 0 1 3 5) (66 27 1) simpleGrading (1 1 1)
    hex (4 2 6 8 5 3 7 9) (66 21 1) simpleGrading (1 1 1)
    hex (2 10 12 6 3 11 13 7) (90 21 1) simpleGrading (1 1 1)
    hex (14 15 16 10 14 15 17 11) (90 27 1) simpleGrading (1 1 1)
    hex (10 16 18 12 11 17 19 13) (90 21 1) simpleGrading (1 1 1)
);

edges
(
);

boundary
(
    inlet
    {
        type patch;
        faces
        (
            (0 4 5 0)
            (4 8 9 5)
        );
    }
    outlet
    {
        type patch;
        faces
        (
            (17 19 18 16)
            (15 16 17 15)
        );
    }
    bottom
    {
        type empty;
        faces
        (
            (0 1 1 0)
            (14 15 15 14)
        );
    }
    top
    {
        type patch;
        faces
        (
            (8 6 7 9)
            (6 12 13 7)
            (12 18 19 13)
        );
    }
    obstacle
    {
        type wall;
        faces
        (
            (1 2 3 1)
            (2 3 11 10)
            (14 10 11 14)
        );
    }
    wedge1
    {
        type wedge;
        faces
        (
            (0 1 2 4)
            (4 2 6 8)
            (2 10 12 6)
            (10 16 18 12)
            (14 15 16 10)
        );
    }
    wedge2
    {
        type wedge;
        faces
        (
            (0 5 3 1)
            (5 9 7 3)
            (3 7 13 11)
            (11 13 19 17)
            (14 11 17 15)
        );
    }
);

mergePatchPairs
(
);

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


This case now works OK. Although there is some strange upwind-moving shock wave (caused by imperfect boundaries, I believe). So the problem is probably somewhere in gmsh grid. For further testing, I'll try to generate uniform grid in gmsh and test the case with it.

[OlegSutyrin]

Another update:
Exactly the same mesh as in blockMeshDict, but generated in gmsh using transfinite lines (no surface at r=0.001 now) = solving OK.
Nearly the same gmsh mesh (but without transfinite lines) with almost uniform cell size = solving OK.
New non-uniform gmsh mesh = solving OK.

In conclusion, I think there are two possible answers here:
1. Problem may have been in my naive cutting out of the symmetry axis (i. e. using r=0.001 as bottom edge). It generates very thin wall-boundary faces that successfully pass checkMesh utility, but still lead to some instabilities.
2. Or the problem was in some mesh irregularity near the outlet (see the bottom right corner of step2Cyl.jpg) - new non-uniform mesh was generated slightly different and hadn't this irregularity. To my intuition, such irregularity shouldn't play any role in the area that is at rest (U =(0 0 0)), but who knows...

After this testing I have some ideas on how to approach mesh generating and initial conditions for my more complex simulations. So posting this here wasn't useless after all, thank you guys