CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Programming & Development (https://www.cfd-online.com/Forums/openfoam-programming-development/)
-   -   The problem of changing k value in the first grid near wall. (https://www.cfd-online.com/Forums/openfoam-programming-development/226308-problem-changing-k-value-first-grid-near-wall.html)

qi.yang@polimi.it April 23, 2020 10:18

The problem of changing k value in the first grid near wall.
 
Hi guys,


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

kc=(U\tau)^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::DimensionedField <double, Foam::volMesh>::<anonymous>.Foam::Field<double>::< anonymous>.Foam::List<double>::<anonymous>.Foam::U List<double>::operator[](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 April 23, 2020 13:55

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

qi.yang@polimi.it April 23, 2020 14:40

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 April 23, 2020 16:11

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

qi.yang@polimi.it April 23, 2020 16:50

4 Attachment(s)
I attached my wall functions:(1)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 April 23, 2020 17:40

You've got a copy of the field here
Code:

volScalarField k = tk();
But you actually need a reference to it
Code:

volScalarField& k = tk();

qi.yang@polimi.it April 24, 2020 03:14

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 April 24, 2020 05:14

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 April 24, 2020 07:30

Code:

volScalarField& k = tk.ref();

qi.yang@polimi.it April 24, 2020 07:58

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 April 24, 2020 08:01

Quote:

Originally Posted by Zeppo (Post 767167)
Code:

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 April 24, 2020 13:49

Code:

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


qi.yang@polimi.it April 26, 2020 11:52

Quote:

Originally Posted by Zeppo (Post 767272)
Code:

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 April 26, 2020 18:01

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

Tibo99 June 15, 2020 17:06

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,


All times are GMT -4. The time now is 01:46.