[egerland]

Hello,

I try to simulate the flow around a cube in a 3D-channel with simpleFoam for incompessible flow, steady state.

set up:
inlet and outlet: cyclic
front and back of the channel: symmetryPlane
top and bottom of the channel: wall

with increasing iterations the streamwise velocity decreases. because of this the massflow is decreasing, too, I think (inkompressible steady state flow)?!?

How can I implement a constant mass flow in simpleFoam.

Thanks

Benni

[ata]

Hi Benni
Your problem is ill conditioned mathematically (has no unique solution in steady state). I think you can use a velocity inlet boundary condition and a value of pressure at outlet.
Good luck
Best regards

Ata

[ata]

Hi Benni
Your problem is ill conditioned mathematically (has no unique solution in steady state). I think you can use a velocity inlet boundary condition and a value of pressure at outlet.
Good luck
Best regards

Ata

[egerland]

Hi Ata,
thanks for your answer. But I need cyclic inlet outlet boundary conditions, so I can simulate an array of cubes. For only one cube I tried velocity inlet BC and a pressure BC at the outlet. It works better.
But I`m interested in a constant mass flow for simple foam?!
Thanks
Benni

[herbert]

Hi Benni,

I think what Ata wanted to say is that you must not use cyclic boundary conditions in your case. Just imagine what will happen to your pressure. When you force it to be equal on inlet and outlet, you suppress pressure loss. But without pressure loss, there can be no flow!

Two ways to solve it:

1. Modify your solve to handle pressure gradient explicitly to enable the use of cyclic boundaries (already discussed here)
2. Use directMapped BC for velocity (you can set a fixed average value here to fix mass flux) and for turbulence properties at inlet and zeroGradient at outlet, but use zeroGradient for pressure on both sides. It can be solved with unmodified simpleFoam.

I've tested both ways, producing nearly the same results for a similar problem.

Regards,
Stefan

[egerland]

Hi Stefan,

thanks a lot. now I unterstand the problem. I will try it directly
best regards
Benni

[egerland]

Hi Stefan,

I tried to modify my solver: I found more information here: http://www.cfd-online.com/Forums/ope...implefoam.html

UEqn.H

// Solve the Momentum equation

tmp<fvVectorMatrix> UEqn
(
fvm::div(phi, U)
+ turbulence->divDevReff(U)
==
flowDirection*gradP
);

UEqn().relax();

eqnResidual = solve
(
UEqn() == -fvc::grad(p)
).initialResidual();

maxResidual = max(eqnResidual, maxResidual);


pEqn.H:

p.boundaryField().updateCoeffs();

volScalarField AU = UEqn().A();
volScalarField rUA = 1.0/UEqn().A();
U = UEqn().H()/AU;
//UEqn.clear();
phi = fvc::interpolate(U) & mesh.Sf();
adjustPhi(phi, U, p);

// Non-orthogonal pressure corrector loop
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
fvm::laplacian(1.0/AU, p) == fvc::div(phi)
);

pEqn.setReference(pRefCell, pRefValue);
// retain the residual from the first iteration
if (nonOrth == 0)
{
eqnResidual = pEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
}
else
{
pEqn.solve();
}

if (nonOrth == nNonOrthCorr)
{
phi -= pEqn.flux();
}
}

# include "continuityErrs.H"

// Explicitly relax pressure for momentum corrector
p.relax();

// Momentum corrector
U -= fvc::grad(p)/AU;
U.correctBoundaryConditions();

// Correct driving force for a constant mass flow rate

// Extract the velocity in the flow direction
dimensionedScalar magUbarStar =
(flowDirection & U)().weightedAverage(mesh.V());

// Calculate the pressure gradient increment needed to
// adjust the average flow-rate to the correct value
dimensionedScalar gragPplus =
(magUbar - magUbarStar)/rUA.weightedAverage(mesh.V());

//volScalarField rUA = 1.0/UEqn.A();
U += flowDirection*rUA*gragPplus;

gradP += gragPplus;

Info<< "Uncorrected Ubar = " << magUbarStar.value() << tab
<< "pressure gradient = " << gradP.value() << endl;

Its working, but uncorrected Ubar is osszilating around Ubar (defined in transportProperties).
pressure gradient is increasing, too.

so I get no convergence.



And I have some problems with directMapped BC:


boundary:

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

dictionaryReplacement
{
boundary
{
inlet
{
type directMappedPatch;
sampleMode nearestPatchFace;
offset (0 0 0);
samplePatch outlet;
//sampleRegion region0;

//samplePatch none;
}
}
}


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

error message:

ggert@linux33:~/OpenFOAM/eggert-1.7.0/run/directMapped_test2$ simpleFoam
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.7.0 |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : 1.7.0-21131bcbd876
Exec : simpleFoam
Date : Jul 27 2010
Time : 17:04:45
Host : linux33
PID : 25270
Case : /home/eggert/OpenFOAM/eggert-1.7.0/run/directMapped_test2
nProcs : 1
SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0

Reading field p

Reading field U

Reading/calculating face flux field phi

Selecting incompressible transport model Newtonian
Selecting RAS turbulence model kEpsilon
#0 Foam::error:rintStack(Foam::Ostream&) in "/software/OpenFOAM/OpenFOAM-1.7.0/lib/linuxGccDPOpt/libOpenFOAM.so"
#1 Foam::sigFpe::sigFpeHandler(int) in "/software/OpenFOAM/OpenFOAM-1.7.0/lib/linuxGccDPOpt/libOpenFOAM.so"
#2 Uninterpreted:
#3 Foam::directMappedFixedValueFvPatchField<double>:: updateCoeffs() in "/software/OpenFOAM/OpenFOAM-1.7.0/lib/linuxGccDPOpt/libfiniteVolume.so"
#4 Foam::fvPatchField<double>::evaluate(Foam::Pstream ::commsTypes) in "/software/OpenFOAM/OpenFOAM-1.7.0/applications/bin/linuxGccDPOpt/simpleFoam"
#5 Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>::GeometricBoundaryField::evaluate() in "/software/OpenFOAM/OpenFOAM-1.7.0/applications/bin/linuxGccDPOpt/simpleFoam"
#6 Foam::incompressible::RASModels::kEpsilon::kEpsilo n(Foam::GeometricField<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh> const&, Foam::transportModel&) in "/software/OpenFOAM/OpenFOAM-1.7.0/lib/linuxGccDPOpt/libincompressibleRASModels.so"
#7 Foam::incompressible::RASModel::adddictionaryConst ructorToTable<Foam::incompressible::RASModels::kEp silon>::New(Foam::GeometricField<Foam::Vector<doub le>, Foam::fvPatchField, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh> const&, Foam::transportModel&) in "/software/OpenFOAM/OpenFOAM-1.7.0/lib/linuxGccDPOpt/libincompressibleRASModels.so"
#8 Foam::incompressible::RASModel::New(Foam::Geometri cField<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh> const&, Foam::transportModel&) in "/software/OpenFOAM/OpenFOAM-1.7.0/lib/linuxGccDPOpt/libincompressibleRASModels.so"


Any ideas?

Thanks

Benni

[Akuji]

Hi everyone,
I'm trying to simulate a turbulent flow around the cube. I use a simpleFoam, but have a strange results: stream lines look like laminar flow. But it is impossible, because inlet velocity is 15 m/s.
How can it be?

[ata]

Hi Arina
As you know laminar or turbulent flow is charaterized by Re number. So, how much is your Re number?
If you are sure that your case is turbulent please explain more about your case, solver, turbulence model and boundary conditions.
Good luck
Best regards

Ata

[Akuji]

Originally Posted by ata

Hi Arina
As you know laminar or turbulent flow is charaterized by Re number. So, how much is your Re number?
If you are sure that your case is turbulent please explain more about your case, solver, turbulence model and boundary conditions.

Ata

Hi, Ata
My Re is about 5 000. I use solver simpleFoam-ras,
Boundary conditions are:
inlet - fixedvalue of all parameters,
outlet - zerogradient;
box - fixedvalue of velocity and pressure,
and walls of chanel - zerogradient

[ata]

Hi Arina
Would you please upload some pictures from your geometry and mesh and your setup files?
How much are your residuals?
Good luck
Best regards

Ata

[Akuji]

Originally Posted by ata

Hi Arina
Would you please upload some pictures from your geometry and mesh and your setup files?
How much are your residuals?
Good luck
Best regards

Ata

Em... What is setup files?
Mesh I made in gmsh

[ata]

Hi Arina
I mean "0", "system" and "constant" files.
Is your mesh a 3D one? If it is, it seems that in Y direction you have not enough cells.
Would you please upload some pictures from your "strange results"? And how much are your residuals?
Best regards
Good luck

Ata

[Akuji]

Originally Posted by ata

Hi Arina
I mean "0", "system" and "constant" files.
Is your mesh a 3D one? If it is, it seems that in Y direction you have not enough cells.
Would you please upload some pictures from your "strange results"? And how much are your residuals?
Best regards
Good luck

Ata

Hi Ata
My mesh is a 3D. Why do you think that in Y dir. I haven't enought cells?
Result: you can see, that I have a nonseparated flow

Here are 0,system,constant
http://letitbit.net/download/11314.11ecf6215667cd89018cd6f5f/box.zip.html

[herbert]

Hi Benni,

for your first problem you can try two things:

1. Add underrelaxation for gradP just like

   gradP += gragPplus

2. Increase viskosity and decrease it step by step

For the second problem I've to say that I didn't work with 1.7.0 yet, but the problem has to be inside BC's of k or epsilon (because reading of U works).

Regards,
Stefan

[ata]

Hi Arina
I think that cube affects boundary conditions in Y direction if this boundaries are near the cube. So, you can not assign a fixed value to these boundaries.
I can't download your files. Would you please upload them in an other place such as forshared?
In my opinion from your figure it is not clear that flow does not separated.
Good luck
Best regards

Ata

[Akuji]

Hi ata
Thanks a lot, I found problem and solved it

[egerland]

Hi Stefan,

I already added underelaxation for grad p. Have still problems.
how do I have to set the BC for k and epsilon when I use directMapped. if k (or epsilon) is not fixedValue I get the failure message (see above).

Regards

Benni

[ata]

Originally Posted by Akuji

Hi ata
Thanks a lot, I found problem and solved it

Hi Arina
I am very happy. Would you please tell me how it solved?
Best regards
Good luck

Ata

[Akuji]

Originally Posted by ata

Hi Arina
I am very happy. Would you please tell me how it solved?
Best regards
Good luck

Ata

Well, I;m not sure that it is solved right, but I hope
About boundary condition: at inlet I have fixed value, at box and bottom of
canal I have velocity = 0 and fixed value of pressure (as I remember). At another region I put zero gradient.

So, at result I saw a turbulence whirlwind after box. That is what I want to see)