[qi.yang@polimi.it]

Hi guys,


I would like to change the k value in the first grid near wall by the following formulation:

kc=(U)^2/sqrt(Cmu);

I wrote the code, according to epsilonwallfunction, but I met the problem when I compile it.


derivedFvPatchFields/wallFunctions/kqRWallFunctions/kOngWallFunction/kOngWallFunctionFvPatchScalarField.C: In member function 'virtual void Foam::kOngWallFunctionFvPatchScalarField::updateCo effs()':
derivedFvPatchFields/wallFunctions/kqRWallFunctions/kOngWallFunction/kOngWallFunctionFvPatchScalarField.C:218:120: error: assignment of read-only location '(& k)->Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>::<anonymous>.Foam:imensionedField <double, Foam::volMesh>::<anonymous>.Foam::Field<double>::< anonymous>.Foam::List<double>::<anonymous>.Foam::U List<double>:perator[](celli)'
((nutw[facei] + nuw[facei])*magGradUw[facei])*((nutw[facei] + nuw[facei])*magGradUw[facei])/(Cmu25*Cmu2 );

my codes in kwallfunction:
const scalarField& y = turbModel.y()[patchi];
const tmp<volScalarField> tk = turbModel.k();
const volScalarField& k = tk();
const tmp<scalarField> tnuw = turbModel.nu(patchi);
const scalarField& nuw = tnuw();
const tmp<scalarField> tnutw = turbModel.nut(patchi);
const scalarField& nutw = tnutw();
const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
const scalarField magGradUw(mag(Uw.snGrad()));

forAll(nutw,facei)
{
label celli = patch().faceCells()[facei];

scalar yPlus = Cmu25*y[facei]*sqrt(k[celli])/nuw[facei];



if (yPlus > yPlusLam_)
{


k[celli] =

((nutw[facei] + nuw[facei])*magGradUw[facei])*((nutw[facei] + nuw[facei])*magGradUw[facei])/(Cmu25*Cmu25);

}

[Zeppo]

You declared your k field as a const reference and then you tried to modify it which is not allowed.

[qi.yang@polimi.it]

Thanks for your reply. Yes, you are right. After I changed the codes as following:
tmp<volScalarField> tk = turbModel.k();
volScalarField k = tk();
It can be compiled successfully but the results are absolutely wrong. Could you please give me some suggestions? I also refered the epsilonwallfunction that updates the epsilon near wall cell. But there is always problem I met. Do you think it is possible to write something to update the k near wall cell by using the formulation? In kwallfunction or epsilonwallfunction?

[Zeppo]

Can you show the line in your code where turbModel is defined

[qi.yang@polimi.it]

I attached my wall functions1)kwallfunc (2)epsilonwallfunc

I hope you can help me to check them. I hope I can use one of them to reach my goal. I wanna force k in the near wall center equals to Utau^2/sqrt(Cmu), which I defined in the .C file.

[Zeppo]

You've got a copy of the field here

volScalarField k = tk();

But you actually need a reference to it

volScalarField& k = tk();

[qi.yang@polimi.it]

After I changed like you suggested,
tmp<volScalarField> tk = turbModel.k();
volScalarField& k = tk();

There was compiling error:
derivedFvPatchFields/wallFunctions/kqRWallFunctions/kOngWallFunction/kOngWallFunctionFvPatchScalarField.C: In member function 'virtual void Foam::kOngWallFunctionFvPatchScalarField::updateCo effs()':
derivedFvPatchFields/wallFunctions/kqRWallFunctions/kOngWallFunction/kOngWallFunctionFvPatchScalarField.C:189:28: error: binding reference of type 'Foam::volScalarField& {aka Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>& ' to 'const Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>' discards qualifiers
volScalarField& k = tk();

[qi.yang@polimi.it]

Hi, I changed the codes as following, but the result dosn't change, k in the first grid is still not calculated by my formulation. I cannot understand.

tmp<volScalarField> tk = turbModel.k();
volScalarField k = tk;

const tmp<scalarField> tnuw = turbModel.nu(patchi);
const scalarField& nuw = tnuw();
const tmp<scalarField> tnutw = turbModel.nut(patchi);
const scalarField& nutw = tnutw();

const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];

const scalarField magGradUw(mag(Uw.snGrad()));

const scalar Cmu25 = pow025(Cmu_);

scalarField& kw = *this;
forAll(nutw,facei)
{
label celli = patch().faceCells()[facei];

k[celli] =

((nutw[facei] + nuw[facei])*magGradUw[facei])*((nutw[facei] + nuw[facei])*magGradUw[facei])/(Cmu25*Cmu25);
}

[Zeppo]

volScalarField& k = tk.ref();

[qi.yang@polimi.it]

Thanks for your help. I found that even though I changed the value of k near wall cell, after calculated the kepsilon equations, the k here will be changed again. How to keep the k here always calculated by my formula?

[qi.yang@polimi.it]

Originally Posted by Zeppo

volScalarField& k = tk.ref();

--> FOAM FATAL ERROR:
Attempt to acquire non-const reference to const object from a tmp<N4Foam14GeometricFieldIdNS_12fvPatchFieldENS_7 volMeshEEE>

From function T& Foam::tmp<T>::ref() const [with T = Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>]
in file J:/blueCFD-Core-2017/OpenFOAM-5.x/src/OpenFOAM/memory/tmp/tmpI.H at line 187.

FOAM aborting

[Zeppo]

tmp<volScalarField> tk = turbModel.k();
volScalarField& k = const_cast<volScalarField&>(tk());

[qi.yang@polimi.it]

Originally Posted by Zeppo

tmp<volScalarField> tk = turbModel.k();
volScalarField& k = const_cast<volScalarField&>(tk());

Really thanks! Now this code works however k is still not caculated by my formulation. I found that there is one function in kepsilon turbulence model is used to keep the value of epsilon in near wall cells caculated by the formulation in epsilonwallfunction.
epsEqn.ref().boundaryManipulate(epsilon_.boundaryF ieldRef());
I am thinking that I also need to apply this to k?
Additionally, epsilon_.boundaryFieldRef().updateCoeffs() in kepsilon turbulence model is utillized to update the epsilon and G in the near wall cell while there is no k_.boundaryFieldRef().updateCoeffs() here. How k value in the boundary layer is refered if so?

[Zeppo]

I can't help you much from now on. But your considerations on treating k in analogy with epsilon seem quite reasonable though.

[Tibo99]

Hi @yangqi!

Did you find a solution to your problem?

Because I'm using a different k profile at the Inlet, I was attempting something similar to deal with BC at the wall.

As a workaround, I use a new gradient formulation as a BC at the wall instead. Don't know if it good tough.

By doing a gradient to the formulation you try to implement, it become 0. That why is ok to use zeroGradient at the wall if you use this type of k profile at the Inlet.

Regards,