CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Programming & Development (
-   -   Viscous term expansion (

santiagomarquezd March 23, 2011 08:35

Viscous term expansion
Hi all, I'm studying the implementation of the viscous term in NS equations in the general framework used for example in pisoFoam, nevertheless I can't understand some of the things that are done. Viscous term is:


due incompressibility and applying div to the product, it becomes

laplacian(nuEff, U)+ grad(nuEff)*[(grad(U)+grad(U)^T)] (1)

in laminar regime nuEff=nu and in turbulent regime nuEff=nu+nut. Reading the code we have from pisoFoam.C:


00065            fvVectorMatrix UEqn
00066            (
00067                fvm::ddt(U)
00068              + fvm::div(phi, U)
00069              + turbulence->divDevReff(U)
00070            );

so that the viscous term is contained in turbulence->divDevReff(U) calling. Taking, for example, the kEpsilon.C implementation:


00181 tmp<fvVectorMatrix> kEpsilon::divDevReff(volVectorField& U) const
00182 {
00183    return
00184    (
00185      - fvm::laplacian(nuEff(), U)
00186      - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
00187    );
00188 }

the laplacian term is present, but the second term is completely different to the second one in (1). I tried to derive the above formulation without success and my advisor showed me that both expressions [(1) and which is present in divDevRefff] give different results in deed, using a counterexample with a div(U)=0 field.

All ideas are welcome.


piccinini December 13, 2011 21:44

Hello, Santiago.

Perhaps you have it already solved, but here goes my try:

The laminar part of viscous stress tensor is:
\nabla\cdot\tau = \nabla \cdot \nu [ \nabla U + ( \nabla U )^T ]

The turbulent part of viscous stress tensor is
\nabla\cdot\tau_T = \nabla \cdot \nu_T [ \nabla U + ( \nabla U )^T ]

Now, summing both contributions:
\nabla\cdot(\tau+\tau_T) =  \nabla\cdot[(\nu+\nu_T) \nabla U] + \nabla \cdot [ (\nu+\nu_T) (\nabla U )^T]

The last term in OpenFOAM code is not the same, though:
- it is computed the deviatoric of (\nabla U )^T and not the gradient itself. This is discussed in another topic:

The trace of last term is zero because of the explicit formulation using the velocity field from previous time step, that should be divergence free.

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