CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Running, Solving & CFD (
-   -   Where I am wronG (

lam August 6, 2007 12:20

Hi Hrvoje, Sorry to not rep
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
solve(UEqn() == -fvc::grad(p));
volScalarField AU = UEqn().A();
U = UEqn().H()/AU;
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"

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

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 August 7, 2007 12:15

just take a little look please
just take a little look please

hjasak August 7, 2007 13:06

? You have messed up the signs
? 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.


lam August 8, 2007 03:56

Hi Hrvoje, Thanks to reply
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.



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