51 lines
1.5 KiB
C
51 lines
1.5 KiB
C
|
if (mesh.moving())
|
||
|
{
|
||
|
// The ddt term constructed by hand because it must be Euler
|
||
|
const scalar dt = runTime.deltaT().value();
|
||
|
|
||
|
const scalarField& V = mesh.V().field();
|
||
|
const scalarField& V0 = mesh.V0().field();
|
||
|
scalar sumV = gSum(V);
|
||
|
|
||
|
scalarField volChange = (1.0 - V0/V)/dt;
|
||
|
|
||
|
scalarField divFlux =
|
||
|
- fvc::div(mesh.phi())().internalField();
|
||
|
|
||
|
scalarField conserve = volChange + divFlux;
|
||
|
|
||
|
scalar sumLocalContErr = dt*gSum(mag(conserve*mesh.V()))/sumV;
|
||
|
scalar globalContErr = dt*gSum(conserve*mesh.V())/sumV;
|
||
|
|
||
|
label nMotionErrors = 0;
|
||
|
scalar maxMotionError = 0;
|
||
|
|
||
|
forAll (conserve, cellI)
|
||
|
{
|
||
|
maxMotionError = Foam::max(maxMotionError, mag(conserve[cellI]));
|
||
|
|
||
|
if (mag(conserve[cellI]) > 1e-8)
|
||
|
{
|
||
|
Info<< "Motion conservation error in cell " << cellI << ": "
|
||
|
<< conserve[cellI]
|
||
|
<< " V: " << V[cellI] << " V0: " << V0[cellI]
|
||
|
<< " volume change: " << volChange[cellI]
|
||
|
<< " div flux: " << divFlux[cellI] << endl;
|
||
|
|
||
|
nMotionErrors++;
|
||
|
}
|
||
|
}
|
||
|
|
||
|
if (nMotionErrors > 0)
|
||
|
{
|
||
|
Info<< "Number of motion errors: " << nMotionErrors
|
||
|
<< " out of nCells = " << mesh.nCells() << endl;
|
||
|
}
|
||
|
|
||
|
Info<< "volume continuity errors : "
|
||
|
<< "volume = " << sumV
|
||
|
<< ", max error = " << maxMotionError
|
||
|
<< ", sum local = " << sumLocalContErr
|
||
|
<< ", global = " << globalContErr << endl;
|
||
|
}
|