[mirsaeedi]

Hi guys
I can not register in workshop of openfoam and i am new with openfoam and c++.
My MSC thesis is about solving Supercavitation phenomena with LES method.

is openFoam an apropriate tool to do this thesis?if it is, can i learn openfoam and it's programming on myself. ofcourse i have read the Userguide and programmer's Guide but i think it is not suffiecient for me and my thesis.

thanks
Morteza

[maka]

half-channel simulation.

(1) I ran a full channel computation with channelOodles and wall function. I got good approximate results. The reason that I insist on using a wall function that later I will introduce roughness effects into it.

(2) I'm interested in a flat-plate like boundary layer. So, I ran a case with half-channel and slip boundary condition on topWall. The results were quite different from the full channel (1) using the same setup parameters.

(3) I suspected that maybe, I should not constrain the normal-to-topWall component of U to be zero (as the case with symmetryPlane) and change the b.c. to zeroGrad on the topWall.

(4) The last modification caused the computation to diverge after 50 s but also, continuity error was bigger than usual after the first time step.

Time = 0.2

Mean and max Courant Numbers = 0.0735003 0.350259
BICCG: Solving for k, Initial residual = 0.0506639, Final residual = 6.75043e-06, No Iterations 3
BICCG: Solving for Ux, Initial residual = 0.04048, Final residual = 3.19361e-06, No Iterations 3
BICCG: Solving for Uy, Initial residual = 0.0561796, Final residual = 5.33075e-06, No Iterations 3
BICCG: Solving for Uz, Initial residual = 0.0554753, Final residual = 4.32981e-06, No Iterations 3
ICCG: Solving for p, Initial residual = 0.864096, Final residual = 7.80343e-07, No Iterations 81
time step continuity errors : sum local = 4.30449e-05, global = 4.30371e-05, cumulative = 4.30371e-05
ICCG: Solving for p, Initial residual = 0.048062, Final residual = 9.10684e-07, No Iterations 70
time step continuity errors : sum local = 4.42835e-05, global = 4.42676e-05, cumulative = 8.73047e-05
Uncorrected Ubar = 0.126888 pressure gradient = 0.0335963
ExecutionTime = 1.86 s ClockTime = 1 s

(5) I reduced the time step to 0.1 and the code diverged at the same time (as in (4)).

can any body help me spot the problem? Thanks.

best regards,
Maka

[eugene]

I'm guessing your pressure boundary on the top wall is also zeroGradient. Having both velocity and pressure zeroGradient on the same boundary is not well specified.

Try using a symmetry plane instead

[maka]

your guess is right. Thanks.

Best regards,
Maka

[marhamat]

Hello everybody

I imposed the perturbU package that suggested in this page to the case similar to channel395.
But after runnig perturbation die out fast and i didn't see turbulent structures .
Do anybody have experience in this field

Best Regards
Marhamat

[maka]

channel395 1-2D

I'm trying to run channelOodles in laminar mode on 1 or 2 D geometry. I use this as a quick test when I modify the solver.

If we take the tutorial case channel395 and change no of cell in x and z direction to 3 the case gives the following error:

Time = 0.4

Mean and max Courant Numbers = 0.000384804 0.020025
Foam::error::printStack(Foam:stream&)
Foam::sigSegv::sigSegvHandler(int)
/lib64/tls/libc.so.6 [0x3072c2e2b0]
Foam::BICCG::solve()
Foam::fvMatrix<double>::solve(Foam::Istream&)
Foam::lduMatrix::solverPerformance Foam::solve<double>(Foam::tmp<foam::fvmatrix<doubl e> > const&)
Foam::LESmodels::oneEqEddy::correct(Foam::tmp<foam ::geometricfield<foam::tensor< double>, Foam::fvPatchField, Foam::volMesh> > const&)
Foam::LESmodel::correct()
channelOodles [0x417449]
__libc_start_main
__gxx_personality_v0
Segmentation fault

increasing the number of cells to 10 in x and z, make the solver work with no probelm.

the problem is that to reproduce a 1D problem one needs to use 3 cell in x and z direction with slip b.c. My question, what does that error mean?

I'm aware of the empty batch option Message.

I suspected that it is related to aspect ratio of the cell. To check that I reproduced the same aspect ratio of the 10-cell case, that worked, in a 3-cell case but again, I got an error.

best regards,
Maka

[hjasak]

This looks pretty nasty: are you sure the code is doing what you expect it to do. Specifically, your error is reported in:

Foam::LESmodels::oneEqEddy::correct

and you said you are using a laminar model. How come? (this should help you work it out - it is not meant to be an ironic question)

Hrv

[maka]

This error message was generated using oneEqEddy as in the standard tutorial case. Also using cyclic b.c. The only change was no. of cell in x,z = 3.

If I make the LESmodel=laminar the code diverge. If I make the time step smaller, 0.02, I get the following error:

Time = 0.16

Mean and max Courant Numbers = 3.85052e-05 0.00200338
BICCG: Solving for Ux, Initial residual = 0.0265804, Final residual = 5.12846e-09, No Iterations 2
BICCG: Solving for Uy, Initial residual = 0.762355, Final residual = 1.26971e-07, No Iterations 2
BICCG: Solving for Uz, Initial residual = 0.90839, Final residual = 2.03159e-07, No Iterations 2
ICCG: Solving for p, Initial residual = 0.42383, Final residual = 9.08821e-07, No Iterations 10
*** glibc detected *** corrupted double-linked list: 0x00000000006db9e0 ***
Aborted

I posted the oneEqEddy error since It contained more details. Thanks for your prompt reply.

best regards,
Maka.

[hjasak]

Just to get things clear: all you did is to use the standard channelOodles tutorial with the number of cells in x and z direction set to 3, as follows:

hex (0 1 3 2 6 7 9 8) (3 25 3) simpleGrading (1 10.7028 1)
hex (2 3 5 4 8 9 11 10) (3 25 3) simpleGrading (1 0.0984 1)

If I do that, I get the following message:

--> FOAM FATAL ERROR : index 1001 out of range 0 ... 449

From function
UList<t>::checkIndex(const label)
in file /home/hjasak/OpenFOAM/OpenFOAM-1.3/src/OpenFOAM/lnInclude/UListI.H at line 107.

Chasing it (a bit) further, it says:

Program received signal SIGABRT, Aborted.
[Switching to Thread 1093481520 (LWP 17231)]
0xffffe410 in __kernel_vsyscall ()
(gdb) where
#0 0xffffe410 in __kernel_vsyscall ()
#1 0x40f877d0 in raise () from /lib/libc.so.6
#2 0x40f88ea3 in abort () from /lib/libc.so.6
#3 0x40bd54b0 in Foam::error::abort (this=0x40dfcc00) at error.C:460
#4 0x0807cf46 in Foam::operator<<> (os=@0x8110630, m={fPtr_ = {__pfn = 0x40bd51f4 <foam::error::abort()>, __delta = 0}, err_ = @0x40dfcc00}) at errorManip.H:85
#5 0x0808ad20 in Foam::UList<double>::checkIndex (this=0xbfbbea28, i=1001) at UListI.H:107
#6 0x0808b34e in Foam::UList<double>::operator[] (this=0xbfbbea28, i=1001) at UListI.H:121
#7 0x080c20cd in Foam::fvMatrix<double>::setReference (this=0xbfbbe9cc, cell=1001, value=@0xbfbbf1d0) at fvMatrix.C:509
#8 0x0807a1b8 in main (argc=Cannot access memory at address 0x434f
) at channelOodles.C:102

This is because the reference cell for the pressure in the tutorial is set to 1001 and your small mesh has got only 450 cells. When I change the reference pressure location to fit the actual mesh used, all works fine.

Can you reproduce this or is it one of my bug fixes that sorts things out?

Hrv

[maka]

If I do the following:

1) modify the blockMeshDict as you did.
2) rm -r 0; cp -r 0.org 0;
3) run

I get the following:

Time = 2.8

Mean and max Courant Numbers = 11603.3 241614
BICCG: Solving for k, Initial residual = 1, Final residual = 9.01383e-06, No Iterations 30
BICCG: Solving for Ux, Initial residual = 0.997261, Final residual = 6.32563e-06, No Iterations 13
BICCG: Solving for Uy, Initial residual = 0.916935, Final residual = 3.92525e-06, No Iterations 12
BICCG: Solving for Uz, Initial residual = 0.996409, Final residual = 6.47291e-06, No Iterations 13
ICCG: Solving for p, Initial residual = 0.999871, Final residual = 8.5682e-07, No Iterations 10
time step continuity errors : sum local = 0.0125487, global = -4.92828e-14, cumulative = -1.33652e-05
ICCG: Solving for p, Initial residual = 0.261537, Final residual = 3.32398e-07, No Iterations 10
time step continuity errors : sum local = 0.00860767, global = -5.98098e-15, cumulative = -1.33652e-05
Uncorrected Ubar = -24.9719 pressure gradient = 2.22217e+07
Foam::error::printStack(Foam:stream&)
Foam::sigSegv::sigSegvHandler(int)
/lib64/tls/libc.so.6 [0x3072c2e2b0]
Foam::cyclicFvPatchField<double>::cyclicFvPatchFie ld(Foam::cyclicFvPatchField<do uble> const&, Foam::Field<double> const&)
Foam::cyclicFvPatchField<double>::clone(Foam::Fiel d<double> const&) const
Foam::GeometricField<double,>::GeometricBoundaryFi eld::GeometricBoundaryField(Fo am::Field<double> const&, Foam::GeometricField<double,>::GeometricBoundaryFi eld const&)
Foam::GeometricField<double,>::GeometricField(Foam ::GeometricField<double,> const&)
Foam::GeometricField<double,>::GeometricField(Foam ::GeometricField<double,> const&)
Foam::GeometricField<double,>::GeometricField(Foam ::GeometricField<double,> const&)
Foam::GeometricField<double,>::New(Foam::IOobject const&, Foam::tmp<foam::geometricfield<double,> > const&, Foam::dimensionSet const&)
Foam::tmp<foam::geometricfield<double,> > Foam::sqrt<foam::fvpatchfield,>(Foam::tmp<foam::ge ometricfield<double,> > const&)
Foam::LESmodels::GenEddyVisc::epsilon() const
channelOodles [0x4189eb]
__libc_start_main
__gxx_personality_v0
Segmentation fault

I think you have something different that allowed for such a clear error message:

"--> FOAM FATAL ERROR : index 1001 out of range 0 ... 449"

If I change the refCell, everything works OK. Thanks for your help.

best regards,
Maka.

[hjasak]

No, I am just using the debug version (always do for development work). I think you are suffering from memory violation because you are accessing an array out of bounds.

Now edit channel395/system/fvSolution

and change the reference pressure location under PIS):

pRefCell 1;

and try again.

Eagerly awaiting :-)

Hrv

[maka]

Yes, it is a pRefCell problem. Now, it works. Sorry I do not know how to work in debug mode. Thanks for your help.

Best regards,
Marwan.

[marhamat]

Hi

I'm interested in LES of turbulent flat plate
Do you have any experience in this field ?
What is the challenge in this case?
With which conditions you get right results?

Best Regards
Marhamat

[maka]

Hi Marhamat,

Thanks for the email, that you sent me to my personal email. It is better to conduct the discussion here in the forum, so you can get comments from different people and such comments remain there to help other new users. I will attempt to answer the following quoted message from you below.

I am interested in turbulent flat plate simulation.
So i have some questions in this field.
I hope that i can implement your experience in this field.
Q1)What is the challenge in this case?
Q2)What was your conditions for getting right results?
- How much is the Re and the grid size?
-How long is the case run time?
- How you initialize this case?
- what is your reference for result validating?
And any recommend that you think is impotant for this problem is pleasent for me.

First, as you can see from my messages, I'm a beginner myself but I will try to answer.

Q1) that is too general question, because the answer depends on your background. So, if the case is important for you, start tackling it and try to get help on the problems that you face along your way.

Q2)

- Re number: That is decided by the application for which you want to do such flat-plat computation. Grid size depends on Re number.

- How long: this depends on the physical processes you are modeling so, if you have heat transfer, it may be different than isothermal flow. If you are interested in fully developed flow, then you can check convergence of the statistics by sampling after you reach fully developed state long enough till the statistics aproximatly does not change.

- i.c.: If you are interested in fully developed flow, you may search the forum for "perturbU", which is a utility that may help you to initialize the flow. You may need to modify to match your case.

- validate: look on the web for Direct Numerical Simulation (DNS) result that has Re number high enough to be close to your objective. If you are doing very high Re number computation, you will not find exact match for your case. But it may be useful to try to reproduce similar results to the DNS data by trying to run at the same Re number of it first, as an example.

It may be useful to start studying channelOodles tutorial as a step towards doing flat plate. There is much to be learnt from that tutorial, that will help you in your case. Also it is a good practice that if you find a problem in your specific case, to try to see if the problem is also reproduced in the tutorial case. This make it much easier for other people to help you, without having to do so much effort. Good luck!

Best regards,
Maka.

[marhamat]

Hi Maka
Thanks alot for your kindness

As you see in this page i had some experience in using oodles and channelOodles tutorial for turbulent pipe flow simulation.
The process for getting right result isn't simple&clear.
So i want to know did you get right result in this case?
If OK , please explain the condition of your case and the process for getting right result with more details .

Best Regards
Marhamat

[maka]

The tutorial case of channelOodles gives a good result. You may start there.

Best regards,
Maka

[marhamat]

Hi Maka
I used channelOodles for this case(turbulent flat plate) with below settings.
Dimentions (1000 100 1)mm for 2D analysis with two blocks
(
hex (0 1 2 3 4 5 6 7) (20 100 1) simpleGrading (1 1000 1)
hex (1 8 9 2 5 11 10 6) (80 100 1) simpleGrading (1 1000 1)
);
In block one bed is symmetryplane

Inlet Velocitt(U)=10m/s & nu=1e-6
You see distribution of U after runnig

I don't see any turbulence structures in my case.
What is your idea about this?

Best Regards
Marhamat

[cedric_duprat]

Hi everyone

I would like there to continue the upper discution between Eugene de Villiers and Rolando Maier (April 2005...)
I'm interested in implementation of the pres-grad-sens wall function.don't worry, I want to implement it by myself. :-)
So, I start to read the "nuSgsWallFunction" which is the Spalding model law of the wall (1961).
I'm quite suprised that the wall fonction is independant of y+. I think that following the value of y+, we can be able to know where the first grid point is and then which sort of law of the wall we will use for example:
Werner Wengle:
y+<11.8> u+ = y+
y+>11.8 ===> u+ = a (y+)**(1/7)

so my questions are:

+ what is the goal of the loop : if(utau > VSMALL)
and what is exactly VSMALL.
is it fonction of the first grid points?

+ in the scalar f definition (in "nuSgsWallFunction"), and with Spalding law of the wall :

y+ = u+ + (1/E) * [ exp(K*u+) - 1 - K*u+ -(K*u+)**2/2 - (Ku+)**3/6 - (K*u+)**4/24 ]

scalar f= (utau/(ry[facei]*nuw[facei])) + y+
what is the term (utau/(ry[facei]*nuw[facei]))

+ Can we use the "nuSgsWallFunction" to implement a wall fonction for LES model ? or is there any easyier way ?

+ what about adding a flux (tau_wall) a the first node wheras modifying the viscosity on the wall ? is it the same thing ?

thank you for helping,

Cedric

[eugene]

How did you figure out the wall model is independent of y+? This is certainly not the case.

The nuSgsWallFunction uses Newton-Raphson iteration to extract utau from U and y - i.e. equilibrium wall function. VSMALL is a very small value 1e-30.

If you want to implement a flux boundary condition wall function, you will have to modify all the top level codes and/or turbulence models that use it. Currently wall functions are bundled with viscosity, because it makes calculating the div(tau) much easier.

Just a note, I did try to implement pressure sensitised wall functions for LES during my thesis and it idn't work. The effect was negligible. Of course I might have made a mistake with my implementation. I suggest you use pMean to calculate the pressure gradient instead of instantaneous p.

[eugene]

Oh yes,

Tauw = sqr(utau)/rho

y+ = y*utau/nu
u+ = u/utau

If you use u+ = a (y+)**(1/7), where do you get utau from before you have calculated the wall shear? I simply iterate until I find a converged value for utau.