[turnow]

Hi all,

I am currently trying to implement a Model called "Approximate deconvolution model (ADM)" for wall bounded flows. The ADM model for large eddy simulations has been proposed by S. Stolz and N.Adams. The main ingredient is an approximation of the nonfiltered Field by truncated series expansion of the inverse Filter operator.

Therefore the momentum equation in the channeloodles solver becomes explicit.
After running the code for 10 - 15 timesteps (deltaT = 0.5e-03) the simulation becomes instable and i get pressure peaks in the corners of the discretized channel.

Here is modified channelOodles solver
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //


Info<< "\nStarting time loop\n" << endl;

for(runTime++; !runTime.end(); runTime++)
{
Info<< "Time = " << runTime.timeName() << nl << endl;


# include "CourantNo.H"

# include "readPISOControls.H"

//sgsModel->correct();

autoPtr<lesfilter> gridFilterPtr_
(
LESfilter::New(U.mesh(),sgsModel->LESmodelProperties())
);

LESfilter& gridFilter_(gridFilterPtr_());
dimensionedScalar Xsi
(
"Xsi",
dimensionSet(0, 0, -1, 0, 0),
1.0
);

//the deconvolved field

Ustar = 3.0*(U) - 3.0*(gridFilter_(U)) + gridFilter_(gridFilter_(U));

fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvc::div(gridFilter_(Ustar * Ustar)(), "div(phi,U)")
- fvc::div(nu*dev(symm(fvc::grad(U))), "div(diss,U)")
- Xsi*(gridFilter_(Ustar) - U)
==
flowDirection*gradP
);

if (momentumPredictor)
{
solve(UEqn == -fvc::grad(p));
}

# include "continuityErrs.H"


// --- PISO

volScalarField rUA = 1.0/UEqn.A();

for (int corr=0; corr<nCorr; corr++)
{
U = rUA*UEqn.H();
phi = (fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, U, phi);

for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
fvm::laplacian(rUA,p) == fvc::div(phi)
);

pEqn.setReference(pRefCell, pRefValue);
pEqn.solve();

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

# include "continuityErrs.H"

U -= rUA*fvc::grad(p);

U.correctBoundaryConditions();
}


Please don't look at the programstyle, I am still working on it.

Thanks for your help.
J. Turnow

[david_h]

Have you tried adding some implicitness to the momentum equation ? i.e. "fvm" instead of "fvc" for the viscous term and replacing the relaxation term:
- Xsi*(gridFilter_(Ustar) - U)
with
-Xsi*gridFilter_(Ustar) + fvm::Sp(Xsi,U)

[turnow]

Hi

for the visous term
"grad(U)" is not possible for fvm:: ?
Maybe you have an idea to avoid using fvm::grad(U)?

Furthermore I corrected the convective term in the momentum equation to

+ fvc::div(fvc::interpolate(gridFilter_(Ustar * Ustar)) & mesh.Sf())

I also tried a new explicit pressure correction,
but I got still pressure peaks at the corners of the channel.


volVectorField flux = fvc::div(fvc::interpolate(gridFilter_(Ustar * Ustar)) & mesh.Sf());

for (int corr=0; corr<nCorr; corr++)
{
fvScalarMatrix pEqn
(
fvm::laplacian(p) == -fvc::div(flux + Xsi*(U - gridFilter_(Ustar)))
);

pEqn.setReference(pRefCell, pRefValue);
pEqn.solve();
}

Thanks for all your help.
J. Turnow

[eugene]

Are you using a plane channel test case? If so, try changing the side boundaries to walls.

What happens if you make your timestep 10 times smaller?

I can't see anything obviously wrong, although I must admit I don't know the method.

[david_h]

You are quite right about fvm::grad.

I meant something along lines of the following:
fvm::laplacian(nu,U) + fvc::div(nu*dev(fvc::grad(U)) )
which is similar to the code for the divR() function in "turbulenceModels/incompressible/laminar/laminar.C"

Dave

[turnow]

thanks eugene and dave,

Yes I am using the plane channel test case.
I changed the side boundary's to walls.
But nothing happend. Again my pressure values blow up to -300 in the first 2 timesteps.
By reducing the time step to 0.5e-04 it is the same game.

Also using the implicitness for the viscous term
"- fvm::laplacian(nu, U) - fvc::div(nu*dev(fvc::grad(U)().T()))"
the simulations is getting instable.

I am wondering also about the term

pEqn.setReference(pRefCell, pRefValue);

When using the term my pressure is blowing up to 300, but by neglecting the reference term the pressure rises up to -3.47146e+13.
Another point is the velocity field. Here it is not effected by the high pressure values.
In addition I compared the convective terms from the ADM solver and the original channelOodles solver. The values are in the same range.

If you have access the Journal "Physics of Fluids" and you would like to read more about the ADM model, there's a good articel:
"An approximate deconvolution model for large-eddy simulation with application to incompressibe wall bounded flows" ( Volume 13, Number 4, April 2001// page 997-1015)

[eugene]

Hmm, I seem to recall there was an issue with the reference cell being located adjacent to a cyclic boundary. Try moving the reference cell to somewhere in the middle of the domain.

On the other hand, the fact that it happens at all the corners seems to indicate that the origin of the issue is elsewhere.

If the above doesn't work, you will have to start stripping stuff out until you are left with the minimum system that reproduces the problem. The first step would be to rewrite the standard channelOodles explicitly to see if you get the same issue.