[Fransje]

Well, applying exactly those steps on a fresh pitzDailyDirectMapped tutorial taken from the ./OpenFOAM-1.7.x/tutorials folder, I definitely get different results.

May I suggest you re-copy a fresh tutorial from that folder, and re-run the steps we went through together? It should only take you a few minutes, and I'm curious to see whether you get the same results than what you posted earlier.
Because normally, the settings you used should not give you those pictures.

Kind regards,

Francois.

[panda60]

Dear Francois,

Thanks for your answer so quickly. the time tween Netherlands and Japan is different, so yesterday night I went home to sleep.
I see the directMapped utility was modified again last month, it seems that a lot have been changed. So I will update my version again, and then check again.

Thanks.
2010-12-08


BUG: directMappedPatchBase : split off nearInfo class

mattijs (author)
December 08, 2010

[panda60]

Dear Francois,

I updaed to the new version, this also happens. Do you think my setting is right ? what is the meaning of offset ? It is the distance of the inlet plane to
the sampling plane, am I right ?

Thanks.

[Fransje]

Dear Jiang,

Well, I investigated the problem somewhat further, and I came to the conclusion that there is a slight misbehaviour in the directMapped routine.

Normally, to get vortex shedding as you show in your plots, you should have to have the setAverage option put to true.
What this does is that is takes the velocity information from the recycle plane you specify, and rescales it so as to have an average of the value specified by average. This not only makes sure that you get a more or less fixed velocity at the entrance, but it also act as a form of driving force, to "force" your fluid in your domain. If you wouldn't have that "driving force", then the velocity at the inlet would gradually become smaller and smaller, because the velocity at the recycle plane would also become smaller and smaller due to friction and viscous dissipation.

So what researchers normally do when using such a recycle method, is that they have the section between the inlet and the recycle plane driven by a pressure gradient to ensure constant velocity (so the Navier-Stokes equations with a source term), and have the rest of the domain solved normally (N-S equations w/o source term ).

That is why I didn't expect vortex shedding in the examples you showed me. Having setAverge true also means you do not need to use potentialFoam to initialize your flow field.

Now. While re-investigating this test case, I noticed a few things:

1. Using the nearestCell option in ./constant/polyMesh/boundary with setAverage false, the flow actually accelerates in time instead of decelerating.
   This is perhaps due to some interaction between the 2-D vortices (truly unphysical phenomenon) and the flow in the narrow entrance.
2. Using the nearestFace option with setAverage false, the flow dies out immediately.
   I don't really understand why, as I would expect a behaviour somewhat similar to that obtain in the previous bullet.. But I am investigating that right now, for my own work.
3. Using the nearestFace option with directMappedVelocityFixed gives the same result as point 2. Which I also find strange.
4. Cases 1 and 2 behave as expected when setting setAverage false, with an average of 15.

So there might be another unexpected bug in those routines, although a part does work. What are you trying to achieve with these test cases?
Also, please be aware that when using directMappedVelocityFlux to remap information from a plane to the inflow, the routine does not give the correct sign to the phi values, so you will have to change that yourself in the code.

The offset value is the value from the inlet patch to the recycle plane you would like to use.

I also could not access the link you provided. Could you re-post it so I can give it a look?

Kind regards,

Francois.

[panda60]

Dear Francois,
It is indeed interesting. setAverge true may drive the flow. But this average just a uniform value, for example(10,0,0). But a profile average(Ux=f(z)) is what I want. But I don't understand why the flow dies out so immediately if setAverage false was used. It shouldn't like this.

My purpose is to achieve Kataoka method, which is a simplied Lund method, but is more easy to use than Lund method. My case is atmospheric boundary layer simulation, but not duct flow.The coming wind is usually a profile, which is changing with height. And Kataoka's method can achieve this goal. In Kataoka's mehtod, the boundary layer thickness is treated as constant, and a velocity profile which coresponding to the experiment can be insured, only fluctuating part was recycled by mutiplying a damping function. The following is Kataoka's formulation:

To achieve this, do you think I also need pressure gradient to driver the flow?

From tomorrow the students will do the final examination test, anybody was not allowed to go inside the school, so I will in home. But I would like to discuss this next monday.
Thanks.

[panda60]

The link is just the OpenFOAM bugs Reporting web.

http://www.openfoam.com/mantisbt/my_view_page.php

[panda60]

Dear Francois,
It seems that, for "directMappedVelocityFlux ", "nearestFace" only can work if the distance you set is near boundary position. That means, nearestFace only work for boundary patch, but doesn't do interpolation for internal patch. So it nearly the same function like "nearestPatchFace".

Please see bug report.

[panda60]

Dear Francois and James

Please don't stop to discuss here. I need your help. Now it is very clear that,
1. you don't have position to modify "nearestPatchFace" , whereas "nearestFace" already gives you a port to modify. So if you want to achieve recycling purpose, among the three method, "nearestFace" is the only one that can be used.

2. For the method "nearestFace" , it doesn't do interpolation for internal plane, but if it is boundary , everything is OK. I want to know if OpenFOAM has internal patch, then when we generate mesh, we give a name for recycling plane, and then nearestFace can recognize this name, don't need to do interpolation.
But I hope the developer can do this interpolation for us. Because it is difficulty for most users to do this kind of work.

3. I think we can just use the outlet boundary as the recycling plane. Indeed, the internal plane is the best, but really we must use the internal plane ?

[panda60]

Now I try to use outlet boundary as recycling plane. The following is code:
//Kataoka method to generate inflow turbulence
case directMappedPolyPatch::NEARESTFACE:
{
vectorField allUValues(nbrMesh.nFaces(), vector::zero);
//scalarField allPhiValues(nbrMesh.nFaces(), 0.0);
forAll(UField.boundaryField(), patchI)
{
const fvPatchVectorField& Upf = UField.boundaryField()[patchI];
const vectorField& Uapf = UaField.boundaryField()[patchI];
//const scalarField& phipf = phiField.boundaryField()[patchI];
label faceStart = Upf.patch().patch().start();
.............
............
.............
//forAll(Upf, faceI)
//{
// allUValues[faceStart + faceI] = Upf[faceI];
// allPhiValues[faceStart + faceI] = phipf[faceI];
//}

scalar delta = 0.25; // boundary layer thickness
forAll(Upf, faceI)
{
scalar z = faceCentres[faceI].z();
scalar thita = z/delta;
// damping function, which impedes the growth of boundary layer
scalar phithita = 0.5*(1-tanh(8.0*(thita-1.0)/(0.7-0.4*(thita-0.3)))/tanh(8.0));
// predetermined Ux mean velocity profile, usually comes from experimental measurement
scalar Up = -52003*z*z*z*z*z*z+55651*z*z*z*z*z-23271*z*z*z*z+4790.9*z*z*z-509.26*z*z+29.481*z+0.3759;

scalar Ux = Up+phithita*(Upf[faceI].x()-Uapf[faceI].x());
scalar Uy = phithita*(Upf[faceI].y()-Uapf[faceI].y());
scalar Uz = phithita*(Upf[faceI].z()-Uapf[faceI].z());

allUValues[faceStart + faceI] = vector(Ux, Uy, Uz);
}
}
mapDistribute::distribute
(
Pstream::defaultCommsType,
distMap.schedule(),
distMap.constructSize(),
distMap.subMap(),
distMap.constructMap(),
allUValues
);
newUValues.transfer(allUValues);
//mapDistribute::distribute
// (
// Pstream::defaultCommsType,
// distMap.schedule(),
// distMap.constructSize(),
// distMap.subMap(),
// distMap.constructMap(),
// allPhiValues
// );
// newPhiValues.transfer(allPhiValues);
break;
}

The red words part is the most straightforward part. Because this part give us a chance to modify to achieve recycling method.

[panda60]

I give Ux=1.5 as initial field and Inlet patch value. After 7 seconds later, the flow field seems to develop, and have the boundary layer thickness which I want. I will wait some days to have a look at.

In the following, left picture is instantaneous velocity , right is mean velocity. The flow field seems to begain to become fluctuate.

[Fransje]

Dear Jiang,

Sorry for the delay, I was away for family obligations.

Well, apparently the nearestFace option is not applicable, which is a shame. But perhaps you can try modifying the nearestCell variant of the directMappedFixedValue to do what you want on values from the interior of the domain?
Of course, sampling the outlet flow as you do is also a possibility if the geometry is constant throughout. But you have to be aware of the fact that the values you sample at the boundary might be slightly off due to the influence of the boundary conditions, and that it might be better to sample the flow further away from the boundary.

Taking a sample patch from the interior does not work to my knowledge. Or at least, I have tried, but it doesn't work. Foamers, please correct me if I'm wrong.

For readability and to avoid errors, you can perhaps try using the C++ inbuilt pow( ) function instead of having z*z*z*z*z*z. Something like:

scalar Up = -52003*pow(z, 6.0) + 55651*pow(z, 5.0) -23271*pow(z, 4.0) + ...

Since you are using directMappedVelocityFlux as basis for your program, it might also be a good idea to update your fluxes at the end of the programme, by using an interpolation of your new U values. Something like:

phi == ( newU & patch.Sf() );

This will avoid having a mismatch between the new U field you provide, and the "fluxes" which were updated using the previous U field.

This is all what comes to my mind for the time being, as your programme already seems to "work". Or at least, seems to be doing something looking like fluctuations..

Kind regards,

Francois.

[panda60]

Dear Francois,
Thanks for your reply! Indeed, if nearestCell can be used, that will be the best. But I have tried directMappedFixedValue , this boundary file is too complex, and nearly can't be compiled to your own solver, maybe that is why directMappedVelocityFlux file exists, and Engene recommended we use this. But it is a pity that nearestCell doesn't exists in this file.
And I also agree internal plane is the best compared with outlet boundary. But as you know, this work is difficult for user to do.

I think use outlet boundary as recycling plane is OK, I think my flow field already develop, but the kinetic energy is still too small compared with experiment, the following is the figure. I don't know why. Do you know how to make flow field more fluctuating ?

[Fransje]

Dear Jiang,

• What type of grid are you using for you test case?
• What are the first few y+ values of your grid near the wall?
• What is your Reynolds number?
• What is your mean free-flow velocity?
• How long did you let your problem run before making your graph?

Kind regards,

Francois.

[panda60]

Dear Francois,
My mesh is totally rectangular grid, in most regrions, the y+=1. My free-flow velocity is a profile which changes with height. the largest velocity is in above boundary layer height position, is around 1.5 m/s. Simulation time goes on nearly 15 waves. And I check the new result, the kinetic energy nearly doesn't change, I think the turbulence already achieve equilibrium, this maybe the final result.

I think this low fluctuation may come from two reasons:
1. Outlet boundary as recycling plane, and the boundary condition influence the sampling.
2. My damping function is not good , the flow field was damped too much, so the turbulence can't develop too much, I will find help from somebody.

Thanks.

[Fransje]

Dear Jiang,

Well, it looks like your fluctuations are under resolved in the boundary-layer part driven by viscosity. This will then have a major effect on the rest of the kinetic energy in your boundary layer, as the turbulent production will be way too low.

I didn't quite analyse the formula you were implementing, but perhaps you should first try to focus on the damping of the viscous terms.

Kind regards,

Francois.

[janara]

Hi Foam users,

I am using interFoam solver for my two-phase flow simulations in a rectangular channel with low Reynolds RAS model (LaunderSharmaKE).
I have error message during simulation which is given below. Please do let me know how this problem is solved.



-------------------------------------------------------------------
request for volScalarField RASModel::G from objectRegistry region0 failed
available objects of type volScalarField are

17
(
nut
K
((2*nu)*magSqr(grad(sqrt(k))))
rho
k
alpha1_0
nu
(2*magSqr(symm(grad(U))))
nu1
(nut*(2*magSqr(symm(grad(U)))))
rho_0
nu2
p
(((C2*(1-(0.3*exp(-min(sqr((sqr(k)|(nu*epsilon))),50)))))*epsilon)|k)
(((2*nu)*nut)*magSqr(grad(grad(U.component(0)))))
alpha1
epsilon
)
#0 Foam::error:rintStack(Foam::Ostream&) in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/lib/linux64GccDPOpt/libOpenFOAM.so"
#1 Foam::error::abort() in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/lib/linux64GccDPOpt/libOpenFOAM.so"
#2 Foam::Ostream& Foam:perator<< <Foam::error>(Foam::Ostream&, Foam::errorManip<Foam::error>) in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/applications/bin/linux64GccDPOpt/interFoam"
#3 Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const& Foam:bjectRegistry::lookupObject<Foam::Geometric Field<double, Foam::fvPatchField, Foam::volMesh> >(Foam::word const&) const in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/lib/linux64GccDPOpt/libinterfaceProperties.so"
#4 Foam::incompressible::RASModels::epsilonWallFuncti onFvPatchScalarField::updateCoeffs() in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/lib/linux64GccDPOpt/libincompressibleRASModels.so"
#5 Foam::fvMatrix<double>::fvMatrix(Foam::GeometricFi eld<double, Foam::fvPatchField, Foam::volMesh>&, Foam::dimensionSet const&) in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/lib/linux64GccDPOpt/libincompressibleRASModels.so"
#6 Foam::tmp<Foam::fvMatrix<double> > Foam::fvm::Sp<double>(Foam:imensionedField<doubl e, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&) in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/lib/linux64GccDPOpt/libincompressibleRASModels.so"
#7 Foam::incompressible::RASModels::LaunderSharmaKE:: correct() in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/lib/linux64GccDPOpt/libincompressibleRASModels.so"
#8
in "/home/jmarati/OpenFOAM/OpenFOAM-1.6/applications/bin/linux64GccDPOpt/interFoam"
#9 __libc_start_main in "/lib/libc.so.6"
#10
at /usr/src/packages/BUILD/glibc-2.9/csu/../sysdeps/x86_64/elf/start.S:116


From function objectRegistry::lookupObject<Type>(const word&) const
in file /home/dm2/henry/OpenFOAM/OpenFOAM-1.6/src/OpenFOAM/lnInclude/objectRegistryTemplates.C at line 140.

FOAM aborting
-----------------------------------------------------------------------


thanks,
Janara

[maysmech]

Dear Foamers,
I have simulated backward facing step with LES model.
i created geometry with 200,000 grid. it is 3D with cyclic B.C for front and back walls. inlet conditions (U and k) are set same as Driver & Seegmiller paper. outlet pressure is set zero.
I don't understand why:

1- Minimum pressure is -278 m2/s2. as you know it means -278000Pa. it is -2.78 Bar but it is not real to have less than -1 Bar.

2- Reattachment point should be near x=6H which H is step height. but it is more than 11 H in this simulation. what is the problem?

Thanks in advance,

[grandgo]

Originally Posted by Fransje

Dear Yashar,

Yes, you can use DirectMapped method for turbulent flow generation. That's what is was developed for. Whether you can apply it to your jet flow case is an other question.
Maybe you can get away with a random turbulent inlet, as your case is somewhat similar to a flow over a backward facing step, for which the quality of the initial turbulence if of less importance. You should refer to relevant literature over the influence of the quality of turbulent information at the inlet for different flow cases. Think of Lund etal 1998, Keating and Piomelli 2004, etc.

Yes, cyclic BC are used in OpenFOAM. Have a look at the channelFoam tutorial.

Kind regards,

Francois.

hi francois,

you seem to have a lot of experience with directmapped boundary conditions and since i badly need to understand this issue and have no clue, you are my man . i've read the topics in this forum, of course, but it's still not concinving.

FIRST OF ALL:
i have a pipe with diameter D and a length of 5*D. i want a LES and i'm using a modified pisoFoam (energy equation included) with the dynSmagorinsky model. the reynolds number of the pipe ( Re(D) ) is defined, so that the axial velocity has to be fixed.

MY GOAL:

• i want to obtain a fully developed flow through the pipe, either with 'cyclic' or 'directMapped' patchtype. (INLET <---> MAPPING PLANE)

• respectively, i DON'T want to define a flow function at the inlet. the flow should develop on its own.

QUESTIONS:

• the first question that comes to my mind is: if i use a 'mapping plane' in a distance of, lets say, 2*D from the inlet and after a while i note, that the simulation has reached a fully developed flow, how can i set up, that the flow 'leaves the cycle' and just continues with the remaining part of the pipe?
• second question is: do i have to create a patch (in other words: the mapping plane) in the meshing process itself, OR can i define the distance between the directmappedpatch (INLET) and the mapping plane in the boundary files?
• third question: how can i force a constant axial velocity of the flow over the whole length of the pipe, so that the slowdown, resulting from friction, is compensated?

i searched in the forum for satisfying answers (i'm relatively new with OF) but i wasn't successful til now.
so i really hope, that you can help me, francois.

other FOAMers are welcome too of course, to help me...

thanks a lot
best regards
grandgo

[Fransje]

Dear grandgo,

Sorry for my late reply, I've been away for quite a while.
I will try to answer the questions I know something about, but I'm not sure I'll be able to answer them all. Let's see.

• the first question that comes to my mind is: if i use a 'mapping plane' in a distance of, lets say, 2*D from the inlet and after a while i note, that the simulation has reached a fully developed flow, how can i set up, that the flow 'leaves the cycle' and just continues with the remaining part of the pipe?

Well, if I were you, I would try obtaining a fully developed flow using cyclic boundary conditions with a pressure-gradient driven flow. That would be the most accurate way of doing it.
If you want a specific "mean" flow velocity in your domain, you could also try using the directMappedFixedValue boundary condition, with the setAverage flag set to true, in order to rescale your outflow velocity field to an inflow field with a certain maximum velocity. (See the pitzDailyDirectMapped tutorial)

I didn't quite understand what you were trying to do after you obtained a fully developed flow (i.e. let the flow 'leave the cycle') so I'll let that open for now.

• second question is: do i have to create a patch (in other words: the mapping plane) in the meshing process itself, OR can i define the distance between the directmappedpatch (INLET) and the mapping plane in the boundary files?

Well... I'm afraid that the recycling procedures don't allow you to make a copy from a plane within your computational domain, unless you program that yourself... :-P So you'll have to use your outlet plane as recycling plane... It will influence the quality of the flow you are recycling though. It's up you to determine whether the quality is good enough.
The sampling plane will be defined in the boundary file, and should be one of your boundary patches.

• third question: how can i force a constant axial velocity of the flow over the whole length of the pipe, so that the slowdown, resulting from friction, is compensated?
  

I guess adding a (constant) source term to the N.S. should do the trick. Thick of a constant pressure gradient compensating for the friction. You'll have to program that yourself though, but it's not difficult. It's just a question of taking the pisoFoam solver and adding the source term you need to the programmed solver. Something like:

dimensionedScalar gradP = dimensionedScalar("1", dimensionSet(0, 1, -2, 0, 0, 0, 0), 1.0);
vector uDir = vector(1.0, 0, 0);
Info << "gradP-> " << gradP << endl;

fvVectorMatrix UEqn
(
    fvm::ddt(U)
    + fvm::div(phi, U)
    + turbulence->divDevReff(U)
    ==
    uDir*gradP
);

Good luck!

On a final note, be aware of the fact that the dynSmagorinsky model in OpenFOAM is not truly locally dynamic. The coefficients are in fact determined based on domain averages, and not on local conditions. So yes, they are dynamic because they are updated at every time step, but no, they are not changed within your domain to account for local conditions.

Kind regards,

Francois.

[Fransje]

Dear grandgo,

I stand corrected.
You could specify an internal sampling plane using the directMappedFixedValue boundary condition, together with the nearestCell option.

Kind regards,

Francois.