[lam]

Hi Hrvoje,

Sorry to not replying you for the last time (about a_N coefficients. I re-read your thesys, but it stil not ok for me), but in emmergency, I want to know where I'm wrong in the implementation of a Stokes solver:
I begin with the simpleFOam.C file: so I just removed the div(phi,U) and changed some "minus" by "plus", according to he discretisation of Navier-Stokes, explained in your thesys.
(in my simulation, I put turbulence OFF)
****************************
tmp<fvvectormatrix> UEqn
(turbulence->divR(U));
UEqn().relax();
solve(UEqn() == -fvc::grad(p));
p.boundaryField().updateCoeffs();
volScalarField AU = UEqn().A();
U = UEqn().H()/AU;
UEqn.clear();
phi = fvc::interpolate(U) & mesh.Sf();
adjustPhi(phi, U, p);

for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
// I put a "minus" here, because now: Up=[H(U)+grad(p)]/ap, so the discretisation of div(U)=0 gives me the equation of pressure below:
fvm::laplacian(1.0/AU, p) == -fvc::div(phi));
pEqn.setReference(pRefCell, pRefValue); pEqn.solve();
if (nonOrth == nNonOrthCorr)
{
// "+=" instead of "-=" because F = S.Uf = S.[H(U)/ap]f + S.[grad(p)/ap]f
phi += pEqn.flux();
}
}

# include "continuityErrs.H"

p.relax();
// "+=" instead of "-=" because Up=[H(U)+grad(p)]/ap
U += fvc::grad(p)/AU;
U.correctBoundaryConditions();
}

********************
the end of the program is the same.
It seems so easy, but I can't get by myself.
I just need your patience to tell me where I'm wrong.

Thanks in advance,

Lam

[lam]

just take a little look please

[hjasak]

? You have messed up the signs: just delete the convection terms and use a symmetric solver for the momentum equation.

I haven't got time to go through this in detail, but it should really be as simple as that.

Hrv

[lam]

Hi Hrvoje,

Thanks to reply me.
Ok, I'm going to delete the convection term...
even if I don't understand why the signs are wrong.
Hope that you will explain me when you get more time.

Cheers,

Lam