Question about rhoSimpleFoam "if (transonic)"
 

Hi,
I try to study rhoSimpleFoam. In the pEqn.H of rhoSimpleFoam, there is a line "if (transonic)" see below for the file. I understand the function of the that is using incompressible SIMPLE method for low mach number cases and adding convertive term for compressible flow. (by the way, the compressible part in the rhoSimpleFoam of OpenFoam-dev1.5 is different).

I did output the the value transonic by adding
"Info<<"transonic = "<<transonic<<endl;"
in the rhoSimpleFoam.

the output of the is always "0" (see below)

My question are :

1, how to control transonic? it is specified at somewhere or it is defined by mach number?

2, if transonic = 0, it is always solving an incompressible problem but not a compressible problem.

Thank
yu

==============================================

Time = 172
smoothSolver: Solving for Ux, Initial residual = 0.130169, Final residual = 6.88999e-05, No Iterations 2
smoothSolver: Solving for Uy, Initial residual = 0.120254, Final residual = 9.18545e-05, No Iterations 2
DILUPBiCG: Solving for h, Initial residual = 0.0154733, Final residual = 4.72774e-05, No Iterations 1
transonic = 0
GAMG: Solving for p, Initial residual = 0.277657, Final residual = 0.00683687, No Iterations 6
time step continuity errors : sum local = 5.41896, global = 0.830172, cumulative = 3.09577
rho max/min : 1.79102 1.26744
ExecutionTime = 53.24 s ClockTime = 68 s

===============================================





===========================================
below is the modified pEqu file in the folder of rhoSimpleFoam
================================================
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
volScalarField rUA = 1.0/UEqn().A();
U = rUA*UEqn().H();
UEqn.clear();
bool closedVolume = false;
Info<<"transonic = "<<transonic<<endl; // output transonic
if (transonic)
{
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)*(fvc::interpolate(U) & mesh.Sf())
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
fvm::div(phid, p)
- fvm::laplacian(rho*rUA, p)
);
// Relax the pressure equation to ensure diagonal-dominance
pEqn.relax(mesh.relaxationFactor("pEqn"));
pEqn.setReference(pRefCell, pRefValue);
// retain the residual from the first iteration
if (nonOrth == 0)
{
eqnResidual = pEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
}
else
{
pEqn.solve();
}
if (nonOrth == nNonOrthCorr)
{
phi == pEqn.flux();
}
}
}
else
{
phi = fvc::interpolate(rho)*(fvc::interpolate(U) & mesh.Sf());
closedVolume = adjustPhi(phi, U, p);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
fvm::laplacian(rho*rUA, p) == fvc::div(phi)
);
pEqn.setReference(pRefCell, pRefValue);
// Retain the residual from the first iteration
if (nonOrth == 0)
{
eqnResidual = pEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
}
else
{
pEqn.solve();
}
if (nonOrth == nNonOrthCorr)
{
phi -= pEqn.flux();
}
}
}

#include "incompressible/continuityErrs.H"
// Explicitly relax pressure for momentum corrector
p.relax();
U -= rUA*fvc::grad(p);
U.correctBoundaryConditions();
// For closed-volume cases adjust the pressure and density levels
// to obey overall mass continuity
if (closedVolume)
{
p += (initialMass - fvc::domainIntegrate(psi*p))
/fvc::domainIntegrate(psi);
}
rho = thermo.rho();
rho = max(rho, rhoMin);
rho = min(rho, rhoMax);
rho.relax();
Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;

Yu

I also tried to add the option in the fvSolution in the system folder.

SIMPLE
{
nNonOrthogonalCorrectors 0;
pMin pMin [ 1 -1 -2 0 0 0 0 ] 0.100;
rhoMax rhoMax [1 -3 0 0 0 0 0] 100;
rhoMin rhoMin [1 -3 0 0 0 0 0] 1e-3;
transonic true;
}

now output value is transonic is 1; see below for the output;

=================================================
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.6.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : 1.6.x-c9721e2ef326
Exec : rhoSimpleFoamMod
Date : Apr 12 2010
Time : 15:01:50
Host : yu-desktop
PID : 3665
Case : /home/yu/OpenFOAM/yu-1.6.x/run/bumper4/rhoSimpleFoam/vin10
nProcs : 1
SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time
Create mesh for time = 0
Reading thermophysical properties
Selecting thermodynamics package hPsiThermo<pureMixture<sutherlandTransport<specieT hermo<hConstThermo<perfectGas>>>>>
Reading field U
Reading/calculating face flux field phi
Creating turbulence model
Selecting RAS turbulence model laminar
Starting time loop
Time = 1
smoothSolver: Solving for Ux, Initial residual = 1, Final residual = 7.78637e-11, No Iterations 2
smoothSolver: Solving for Uy, Initial residual = 1, Final residual = 4.67024e-07, No Iterations 2
DILUPBiCG: Solving for h, Initial residual = 1, Final residual = 1.53609e-06, No Iterations 1
transonic = 1

--> FOAM FATAL IO ERROR:
Cannot find relaxationFactor for 'pEqn' or a suitable default value.

==============================================



But there is a new error.

Please help.

Yu

Just add 'pEqn .3' entry in your relaxationFactors :)

Still no sure how to use transonic
 

Thanks for your reply.
it can start calculation. but deverged.

Not sure how to use rhoSimpleFoam for a compressible flow.

It is correct to put "transonic" as "true" , like below, or not?

SIMPLE
{
nNonOrthogonalCorrectors 0;
pMin pMin [ 1 -1 -2 0 0 0 0 ] 0.100;
rhoMax rhoMax [1 -3 0 0 0 0 0] 100;
rhoMin rhoMin [1 -3 0 0 0 0 0] 1e-3;
transonic true;
}


In OF-1.5-dev, there is no option for transonic at all. Can some one explain why. Thanks

Yu

Transonic is unstable!
 

'transonic' is a boolean so you can put 'true' or 'yes'.
I have stability problems too with this boolean on :( As transonic equations are strongly non-linear, it seems quite difficult to make a computation converge. I am still investigating, especially trying to understand why density standard deviation (max minus min) is so high, even with both very low CFL number and relaxation factor for rho.
Hope it'll help.:)