CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)
-   -   Conservative form of Navier Stokes equation. (https://www.cfd-online.com/Forums/openfoam-solving/96427-conservative-form-navier-stokes-equation.html)

 balkrishna January 23, 2012 01:43

Conservative form of Navier Stokes equation.

Hi,

I wish to modify the pimpleFoam solver for varying density. As a first step, I decided to change the momentum equation to the conservative form. To do so I changed the momentum equation to
Code:

fvVectorMatrix UEqn
(
fvm::ddt(rho, U) +
fvm::div(phi, U)
- fvm::laplacian(mu, U)
- fvc::div(mu*dev2(fvc::grad(U)().T()))
);

here rho is declared as a volScalarField. The pressure equation is coded in the following way :
Code:

{

Info<<" In the most dreaded pEqn "<<endl;
volScalarField rUA = 1.0/UEqn.A();
surfaceScalarField rhorUAf("(rho*(1|A(U)))", fvc::interpolate(rho*rUA));

U = rUA*UEqn.H();

phi = fvc::interpolate(rho)*
(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rUA, rho, U, phi)
);

surfaceScalarField buoyancyPhi = -rhorUAf*ghf*fvc::snGrad(rho)*mesh.magSf();
phi += buoyancyPhi;

for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix p_rghEqn
(
fvc::ddt(rho) +
fvc::div(phi)
- fvm::laplacian(rhorUAf, p_rgh)
);

p_rghEqn.solve
(
mesh.solver
(
p_rgh.select
(
(
finalIter
&& corr == nCorr-1
&& nonOrth == nNonOrthCorr
)
)
)
);

if (nonOrth == nNonOrthCorr)
{
// Calculate the conservative fluxes
phi += p_rghEqn.flux();

// Explicitly relax pressure for momentum corrector
p_rgh.relax();

// Correct the momentum source with the pressure gradient flux
// calculated from the relaxed pressure
U += rUA*fvc::reconstruct((buoyancyPhi + p_rghEqn.flux())/rhorUAf);
}
}

p = p_rgh + rho*gh;

DpDt = fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);

As a first test I decided to simulate for a case of non changing density for a pipe flow.My rho is specified in the following manner.
Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                |                                                |
| \\      /  F ield        | OpenFOAM: The Open Source CFD Toolbox          |
|  \\    /  O peration    | Version:  1.7.0                                |
|  \\  /    A nd          | Web:      www.OpenFOAM.com                      |
|    \\/    M anipulation  |                                                |
\*---------------------------------------------------------------------------*/
FoamFile
{
version    2.0;
format      ascii;
class      volScalarField;
location    "0";
object      rho;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [1 -3 0 0 0 0 0];

internalField  uniform 1000;

boundaryField
{
outlet
{
type            fixedValue;
value          uniform 1000;
}
inlet
{
type            fixedValue;
value          uniform 1000;
}
wall
{
type            fixedValue;
value          uniform 1000;
}
frontBack
{
type            empty;
}
}

// ************************************************************************* //
0/rho (END)

This gives highly unphysical results as far as the velocity profiles are concerned. Could someone tell me where I might be going wrong in the momentum equation or pressure equation ?

 alberto January 24, 2012 12:22

Take a look at rhoPimpleFoam, and for a solver using p_rgh, you could take a look at buoyantPimpleFoam.

Best,

 balkrishna January 25, 2012 09:33

Solved.

thanks for d help.

 All times are GMT -4. The time now is 19:20.