[Taozi]

4*1 domain, periodic in x and non-slip in y, Ra_cri = 1706, successful with Gerris, unsuccessful with foam-extend-3.0. I saw some people asked about this long ago, but could not find any answer. The problem is there is no instability at all for example at Ra = 1800. I wish to use OpenFOAM for some theoretical analysis, and first of all I have make sure that it can reproduce this classical result. I couldn't figure out what the reason is, maybe grid resolution should be high? maybe run time should be longer? maybe the FVM scheme should be changed? But all seem to be off the point to me.

Another very strange thing is in the following figure (attached, name Selection_003.jpg), it plots the velocity magnitude. I have specified periodic boundary conditions (using cyclic in openfoam) in the x direction, and it doesn't make sense that something (higher velocity magnitude - or just floating point errors) concentrate near the boundary -- there is no horizontal boundary at all!

I have attached the simulation files.

It's four in the morning at my time zone, after fighting against if overnight, I am pretty upset... So could someone offer some help? The case setup is very simple, so reading it and even running it yourself won't take much time, no to mention it's weekend, so please~

Thank you very much!

[wyldckat]

Greetings Tom,

I've finally managed to have a look into this. Problem is that I'm not able to run the simulation case you've provided
I had to create the mesh with blockMesh from OpenFOAM 2.3.x and then used boussinesqBuoyantFoam from foam-extend 3.0 to run the case. But it crashed at 0.2 of the simulation time/iteration.

Therefore, I'm unable to diagnose anything at all
If you have already solved this issue, please provide at least some closure on this topic.
If you have not yet solved this issue, please provide more information and a "working case", in which it's possible to diagnose the same problem you're having.

Best regards,
Bruno

[Taozi]

Hi Bruno, sorry for the late reply.
I also use boussinesqBuoyantFoam from ext-3.0, with slight modification to work with 2.3.x. Attached is the same case file above, together with the solver (inside the case folder), you can wmake (using 2.3.x) in the solver folder (the resultant exe is in the solver folder as well) and
cd ..
./boussinesqBuoyantFoam/boussinesqBuoyantFoam

It runs well and successfully generates the same problem.
Following is the dropbox link:
https://dl.dropboxusercontent.com/u/...eighbenard.zip

Chenguang

[Smed]

Has anyone had a chance to look into this further?

I was able to run Taozi's model from the dropbox attachment using OF 2.3.1. I confirmed that the model looks correct and the boundary conditions should result in a Rayleigh number of 1800, but I don't know why the results don't look right.

[Taozi]

Hi Adam, Thank your for your attention. I still don't have a hint why this happens, especially when I am able to get the correction Ra number with either gerris or fluent. To me this is a serious problem and I had to give up using OpenFOAM for one of my research problem. I will still look into it when I have time, but I wish people with deeper insight into OpenFOAM will come and help :-)

[ssss]

Could you give more informatión about the simulation? Domain, discretization, mesh,etc... Are you sure you are using a correct mesh? You will need to refine a lot to get to the instability to grow if you are near the ra critical

[Smed]

The entire self-contained case is in Taozi's earlier dropbox link.

For a general summary, the 2D cavity is 4 meters wide and 1 meter tall, the mesh is a uniform 128x32, but I tried increasing to 256x128 without any significant effect on results, the temperature difference between plates is 1.8 K, the solution schemes seem to be pretty standard PISO, but again, I'd suggest looking at the system files for more details.

Taozi, can you provide us with your equivalent model for Gerris? I've never used that code, but maybe by looking at the differences we can figure out the problem.

[olivierG]

Hello,

I don't reproduce your result either.
You use transient solver pimple, and don't wait long enough to get the solution.
Here is some pict + case with OF 2.3.x from your case, with same mesh.

regards,
olivier

[Taozi]

Originally Posted by olivierG

Hello,

I don't reproduce your result either.
You use transient solver pimple, and don't wait long enough to get the solution.
Here is some pict + case with OF 2.3.x from your case, with same mesh.

regards,
olivier

Hi, Olivier, thank you very much. I am checking your case. I haven't run it for 2000s yet (still running), but I noticed that the mysterious behavior of U near the side wall is still there (as in the picture of my first post). My understanding is if the domain is periodic and there is no side wall, the solution's behavior should be more or less uniform across the domain. Even instability occurs, it should occur across the domain. The appearance of the solution, however, seems to me that the algorithm some how still feels effect of a virtual wall. That's very strange to me. What do you think of that, did you see it in your simulation results?

[Taozi]

Hi Adam, below is the gerris file I used, you can run it as gerris2D -m filename

Define T1 1.716

4 4 GfsSimulation GfsBox GfsGEdge{}{
        PhysicalParams {L = 1 g = 0}
        Refine 6
        Time{start = 0 end = 1000}

        VariableTracer {} T
        SourceViscosity 0.01
        SourceDiffusion T 0.01
        # beta*g = 0.1
        Source {} V 0.1*T

        # some perturbation is added
        GfsInit {} {
            T = T1 * (0.5 - y) + sin(1000*x)
        }

        # output control
        OutputTime {step=1} stdout
        OutputSimulation {step = 10} ./output-%04.1f.gts
        GfsOutputScalarStats {step = 1} ./maxmin_vel { v = Velocity}
}

# box 0
GfsBox{
    top   = Boundary{
        BcDirichlet T 0
        BcDirichlet U 0
        BcDirichlet V 0
    }
    bottom= Boundary{
        BcDirichlet T T1
        BcDirichlet U 0
        BcDirichlet V 0
    }
}

# box 1
GfsBox{
    top   = Boundary{
        BcDirichlet T 0
        BcDirichlet U 0
        BcDirichlet V 0
    }
    bottom= Boundary{
        BcDirichlet T T1
        BcDirichlet U 0
        BcDirichlet V 0
    }
}
# box 2
GfsBox{
    top   = Boundary{
        BcDirichlet T 0
        BcDirichlet U 0
        BcDirichlet V 0
    }
    bottom= Boundary{
        BcDirichlet T T1
        BcDirichlet U 0
        BcDirichlet V 0
    }
}

# box end
GfsBox{
    top   = Boundary{
        BcDirichlet T 0
        BcDirichlet U 0
        BcDirichlet V 0
    }
    bottom= Boundary{
        BcDirichlet T T1
        BcDirichlet U 0
        BcDirichlet V 0
    }
}
1 2 right
2 3 right
3 4 right
4 1 right
# End of gfs file

[olivierG]

Hello Taozi,

The instability start from the cyclic, but evolve rapidly (solution is already here way before 1000s).
Note that this is without adding initial sinus perturbation like in gerris, so if you want perform a bench test, do the same.

regards,
olivier