[johndeas]

Hi,

I would like to create a random vector field with which I could play with my flow (adding it to velocity components as a white noise, or as a body force in the NS equation,...).

Since I still do not catch well how data is organised within OpenFOAM, could somebody help me ? I know the object for this field would be a volVectorField, but I do not get much...
I would like values to be independent for each positions and components, and spanning between 0 and 1.

Thanks,

JD

[johndeas]

I tried to modify perturbU from De Villier :

#include "fvCFD.H"
#include "Random.H"
#include "wallDist.H"
#include "wallDistReflection.H"
#include "mathematicalConstants.H"

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

int main(int argc, char *argv[])
{
# include "setRootCase.H"

# include "createTime.H"
# include "createMesh.H"

wallDist yw(mesh);
const scalar h = max(yw.internalField());

wallDistReflection reflexVec(mesh);

const volVectorField yDir = reflexVec.n();

IOobject Uheader
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ
);
Info << "Reading U" << endl;
volVectorField U(Uheader, mesh);

IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);

dimensionedVector Ubar
(
transportProperties.lookup("Ubar")
);

vector xDir = Ubar.value()/mag(Ubar.value());

const vectorField& centres(mesh.C());

Random perturbation(1234567); // objet possède la méthode GaussNormal

forAll(centres, celli)
{

scalar devX=1.0 + 0.1*perturbation.GaussNormal();
scalar devY=1.0 + 0.1*perturbation.GaussNormal();
scalar devZ=1.0 + 0.1*perturbation.GaussNormal();

vector zDir = xDir^yDir[celli];
zDir /= mag(zDir);


U[celli] = 3.0*Ubar.value() * (yw[celli]/h - 0.5*sqr(yw[celli]/h));

U[celli] += xDir*devX+yDir*devY+zDir*devZ;

}

U.write();

Info<< endl;

return(0);
}


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

But this line

U[celli] += xDir*devX+yDir*devY+zDir*devZ;

is a source of compilation error:

sandFoam.C:69: error: no match for 'operator+=' in 'U.Foam::GeometricField<foam::vector<double>, Foam::fvPatchField, Foam::volMesh>::<anonymous>.Foam::DimensionedField <foam::vector<double>, Foam::volMesh>::<anonymous>.Foam::Field<foam::vect or<double> >::<anonymous>.Foam::List<foam::vector<double> >::<anonymous>.Foam::UList<t>::operator[] [with T = Foam::Vector<double>](celli) += Foam::operator+(const Foam::tmp<foam::geometricfield<type2,> >&, const Foam::VectorSpace<form,>&) [with Form = Foam::Vector<double>, Cmpt = double, int nCmpt = 3, Type = Foam::Vector<double>, PatchField = Foam::fvPatchField, GeoMesh = Foam::volMesh](((const Foam::VectorSpace<foam::vector<double>, double, 3>&)((const Foam::VectorSpace<foam::vector<double>, double, 3>*)(& Foam::operator*(const Foam::VectorSpace<form,>&, Foam::scalar) [with Form = Foam::Vector<double>, Cmpt = double, int nCmpt = 3](devZ).Foam::Vector<double>::<anonymous>))))'
/home/vinz/OpenFOAM/OpenFOAM-1.4.1/src/OpenFOAM/lnInclude/VectorSpaceI.H:204: note: candidates are: void Foam::VectorSpace<form,>::operator+=(const Foam::VectorSpace<form,>&) [with Form = Foam::Vector<double>, Cmpt = double, int nCmpt = 3]

However I thought I provided my deviation as a three component vector, so I do not understand where the problem is.

[cedric_duprat]

Hi John,

Maybe you can have a look there : src/randomProcesses
there is a file noise and an other turbulence.

I think they could be interresting for you.

With perturbU, you can already add a noise near the wall (without adding the laminar velocity profile for U).

Can you precise what you want do to ?

Regards

Cedric

[johndeas]

I am trying to reproduce results we get on another code (DNS), perturbing a Poiseuille profile with white noise at 1/10 the original local field value.

I noted perturbU and will give it a try later.

Since turbgen generates divergence free turbulent field, I may also try to surimpose one it produces with a Poiseuille profile.

Anyway if someone understands why U[celli] += xDir*devX+yDir*devY+zDir*devZ; throw me an error in the modification I did...

[johndeas]

I corrected my error: yDir was declared of type volVectorField, and, as such, the perturbation had to be written:

U[celli] += devX*xDir+devY*yDir[celli]+devZ*zDir;

instead of:

U[celli] += devX*xDir+devY*yDir+devZ*zDir;