CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM (http://www.cfd-online.com/Forums/openfoam/)
-   -   time step continiuty error (http://www.cfd-online.com/Forums/openfoam/103286-time-step-continiuty-error.html)

MultiphaseFlowsLab June 15, 2012 16:23

time step continiuty error
 
Hi,
I'm using pisoFoam solver and I added energy equation to it. But when I run it at first time step "time stem continuity error" is appeared and in the following time step, my Courant number increases. Do you have any idea why it happens? or how I can find the cell which causes the problem ?

Maryam

niaz June 16, 2012 00:28

Dear Maryam
you should use adjustabletimestep to limit your Co number.
and I suggest you to use pimpleFoam instead of pisoFoam.
it works more stable than pisoFoam.

nimasam June 16, 2012 01:19

some general advices:
1)put your initial time step at 1e-08
2) use " adjustabletimestep to limit your Co number"
if you again encounter problem,
3) check your BC

MultiphaseFlowsLab June 19, 2012 08:57

Thanks guys for your help. It works.
I'm trying to add energy equation to pisoFoam and I followed http://openfoamwiki.net/index.php/Ho...ure_to_icoFoam
When the code is running, it shows me residual for T, but in my results, T is constant.
Do you have any idea?

nimasam June 19, 2012 11:48

whats your BC?

MultiphaseFlowsLab June 19, 2012 12:01

I have two cylinders inside each other, which fluid flows between cylinders. The inner one has constant heat (fixedGradient) and the outer one is zeroGradient. fixedValue for inlet and zeroGradient for outlet.

nimasam June 19, 2012 12:15

it's weird :) does it solve for T at all ?
compare your code with code in heat transfer folder, you will find how to add energy equation

MultiphaseFlowsLab June 19, 2012 12:22

It shows solving, but in any time step, T file shows constant temperature.

nimasam June 19, 2012 12:37

put your code here ;)

MultiphaseFlowsLab June 19, 2012 12:47

Thanks for following up.
I added T field in creatField.H.
and

#include "fvCFD.H"
#include "MysinglePhaseTransportModel.H"
#include "turbulenceModel.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

int main(int argc, char *argv[])
{
#include "setRootCase.H"

#include "createTime.H"
#include "createMesh.H"
#include "createFields.H"
#include "initContinuityErrs.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

Info<< "\nStarting time loop\n" << endl;

while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;

#include "readPISOControls.H"
#include "CourantNo.H"

// Pressure-velocity PISO corrector
{
// Momentum predictor

fvVectorMatrix UEqn
(
fvm::ddt(U)
+ fvm::div(phi, U)
+ turbulence->divDevReff(U)
);

UEqn.relax();

if (momentumPredictor)
{
solve(UEqn == -fvc::grad(p));
}



// --- PISO loop

for (int corr=0; corr<nCorr; corr++)
{
volScalarField rAU(1.0/UEqn.A());

U = rAU*UEqn.H();
phi = (fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, U, phi);

adjustPhi(phi, U, p ,T);

// Non-orthogonal pressure corrector loop
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
// Pressure corrector

fvScalarMatrix pEqn
(
fvm::laplacian(rAU, p) == fvc::div(phi)
);

pEqn.setReference(pRefCell, pRefValue);

if
(
corr == nCorr-1
&& nonOrth == nNonOrthCorr
)
{
pEqn.solve(mesh.solver("pFinal"));
}
else
{
pEqn.solve();
}

if (nonOrth == nNonOrthCorr)
{
phi -= pEqn.flux();
}
}

#include "continuityErrs.H"

U -= rAU*fvc::grad(p);
U.correctBoundaryConditions();


//add these lines...
fvScalarMatrix TEqn
(
fvm::ddt(T)
+ fvm::div(U, T) == fvm::laplacian(DT, T)
);

TEqn.solve();
}
}
//done adding lines...
turbulence->correct();

runTime.write();

Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;

}

Info<< "End\n" << endl;

return 0;
}

nimasam June 19, 2012 13:00

replace it
fvm::ddt(T)
+ fvm::div(U, T) == fvm::laplacian(DT, T)
with
fvm::ddt(T)
+fvm::div(phi, T) == fvm::laplacian(DT, T)

MultiphaseFlowsLab June 19, 2012 13:09

Thanks for checking the code. I have same problem with "phi" also. At first I used "phi" then I changed it to "U":(

nimasam June 19, 2012 13:45

adjustPhi(phi, U, p ,T); ???????

i think this is wrong too! whats the role of p or T ?
and if it was ok ! then change your BC , put two fixed BC for example one at 200 and the other in 500 with no inlet-out, and look whether it works at all!

mm.abdollahzadeh July 10, 2012 11:12

Hi Everybody

Today I started looking at foam again :D.
in pisofoam there is line "Mesh.solver" .

Could any body tell to me what is that?

Best
Mahdi


All times are GMT -4. The time now is 12:23.