ENH: limitVelocity, report count and percent of cells/faces (#2414)

- percent of cells is taken relative to selection size.
- percent of faces is taken relative to the number of boundary faces
  that do not fix velocity themselves.

ENH: avoid correctBoundaryConditions() if values were not limited
This commit is contained in:
Mark Olesen
2022-03-16 15:57:31 +01:00
parent 4495db0302
commit 55f287cd83
4 changed files with 149 additions and 44 deletions

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2021 OpenCFD Ltd.
Copyright (C) 2021-2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -83,41 +83,70 @@ void Foam::fa::limitVelocity::correct(areaVectorField& U)
{
const scalar maxSqrU = sqr(max_);
// Count nTotFaces ourselves
// (maybe only applying on a subset)
label nFacesAbove(0);
const label nTotFaces(returnReduce(faces_.size(), sumOp<label>()));
vectorField& Uif = U.primitiveFieldRef();
for (const label facei : faces_)
{
const scalar magSqrUi = magSqr(Uif[facei]);
auto& Uval = Uif[facei];
const scalar magSqrUi = magSqr(Uval);
if (magSqrUi > maxSqrU)
{
Uif[facei] *= sqrt(maxSqrU/max(magSqrUi, SMALL));
Uval *= sqrt(maxSqrU/max(magSqrUi, SMALL));
++nFacesAbove;
}
}
// Handle boundaries in the case of 'all'
label nEdgesAbove(0);
if (!faceSetOption::useSubMesh())
{
for (faPatchVectorField& Up : U.boundaryFieldRef())
{
if (!Up.fixesValue())
{
forAll(Up, facei)
for (auto& Uval : Up)
{
const scalar magSqrUi = magSqr(Up[facei]);
const scalar magSqrUi = magSqr(Uval);
if (magSqrUi > maxSqrU)
{
Up[facei] *= sqrt(maxSqrU/max(magSqrUi, SMALL));
Uval *= sqrt(maxSqrU/max(magSqrUi, SMALL));
++nEdgesAbove;
}
}
}
}
}
// We've changed internal values so give boundary conditions opportunity
// to correct.
U.correctBoundaryConditions();
// Percent, max 2 decimal places
const auto percent = [](scalar num, label denom) -> scalar
{
return (denom ? 1e-2*round(1e4*num/denom) : 0);
};
reduce(nFacesAbove, sumOp<label>());
reduce(nEdgesAbove, sumOp<label>());
Info<< type() << ' ' << name_ << " Limited "
<< nFacesAbove << " ("
<< percent(nFacesAbove, nTotFaces)
<< "%) of faces, with max limit " << max_ << endl;
if (nFacesAbove || nEdgesAbove)
{
// We've changed internal values so give
// boundary conditions opportunity to correct
U.correctBoundaryConditions();
}
}

View File

@ -62,7 +62,10 @@ void Foam::fv::velocityDampingConstraint::addDamping(fvMatrix<vector>& eqn)
const volVectorField& U = eqn.psi();
scalarField& diag = eqn.diag();
label nDamped = 0;
// Count nTotCells ourselves
// (maybe only applying on a subset)
label nDamped(0);
const label nTotCells(returnReduce(cells_.size(), sumOp<label>()));
for (label celli : cells_)
{
@ -79,10 +82,16 @@ void Foam::fv::velocityDampingConstraint::addDamping(fvMatrix<vector>& eqn)
reduce(nDamped, sumOp<label>());
Info<< type() << " " << name_ << " damped "
// Percent, max 2 decimal places
const auto percent = [](scalar num, label denom) -> scalar
{
return (denom ? 1e-2*round(1e4*num/denom) : 0);
};
Info<< type() << ' ' << name_ << " damped "
<< nDamped << " ("
<< 100*scalar(nDamped)/mesh_.globalData().nTotalCells()
<< "%) of cells" << endl;
<< percent(nDamped, nTotCells)
<< "%) of cells, with max limit " << UMax_ << endl;
}

View File

@ -126,43 +126,53 @@ void Foam::fv::limitTemperature::correct(volScalarField& he)
scalar Tmin0 = min(T);
scalar Tmax0 = max(T);
label nOverTmax = 0;
label nLowerTmin = 0;
// Count nTotCells ourselves
// (maybe only applying on a subset)
label nBelowMin(0);
label nAboveMax(0);
const label nTotCells(returnReduce(cells_.size(), sumOp<label>()));
forAll(cells_, i)
{
const label celli = cells_[i];
if (hec[celli] < heMin[i])
{
nLowerTmin++;
hec[celli] = heMin[i];
++nBelowMin;
}
else if (hec[celli] > heMax[i])
{
nOverTmax++;
hec[celli] = heMax[i];
++nAboveMax;
}
hec[celli]= max(min(hec[celli], heMax[i]), heMin[i]);
}
reduce(nOverTmax, sumOp<label>());
reduce(nLowerTmin, sumOp<label>());
reduce(nBelowMin, sumOp<label>());
reduce(nAboveMax, sumOp<label>());
reduce(Tmin0, minOp<scalar>());
reduce(Tmax0, maxOp<scalar>());
Info<< type() << " " << name_ << " Lower limited "
<< nLowerTmin << " ("
<< 100*scalar(nLowerTmin)/mesh_.globalData().nTotalCells()
<< "%) of cells" << endl;
// Percent, max 2 decimal places
const auto percent = [](scalar num, label denom) -> scalar
{
return (denom ? 1e-2*round(1e4*num/denom) : 0);
};
Info<< type() << " " << name_ << " Upper limited "
<< nOverTmax << " ("
<< 100*scalar(nOverTmax)/mesh_.globalData().nTotalCells()
<< "%) of cells" << endl;
Info<< type() << ' ' << name_ << " Lower limited " << nBelowMin << " ("
<< percent(nBelowMin, nTotCells)
<< "%) of cells, with min limit " << Tmin_ << endl;
Info<< type() << ' ' << name_ << " Upper limited " << nAboveMax << " ("
<< percent(nAboveMax, nTotCells)
<< "%) of cells, with max limit " << Tmax_ << endl;
Info<< type() << ' ' << name_ << " Unlimited Tmin " << Tmin0 << endl;
Info<< type() << ' ' << name_ << " Unlimited Tmax " << Tmax0 << endl;
Info<< type() << " " << name_ << " Unlimited Tmax " << Tmax0 << nl
<< "Unlimited Tmin " << Tmin0 << endl;
// Handle boundaries in the case of 'all'
bool changedValues = (nBelowMin || nAboveMax);
if (!cellSetOption::useSubMesh())
{
volScalarField::Boundary& bf = he.boundaryFieldRef();
@ -183,16 +193,28 @@ void Foam::fv::limitTemperature::correct(volScalarField& he)
forAll(hep, facei)
{
hep[facei] =
max(min(hep[facei], heMaxp[facei]), heMinp[facei]);
if (hep[facei] < heMinp[facei])
{
hep[facei] = heMinp[facei];
changedValues = true;
}
else if (hep[facei] > heMaxp[facei])
{
hep[facei] = heMaxp[facei];
changedValues = true;
}
}
}
}
}
// We've changed internal values so give
// boundary conditions opportunity to correct
he.correctBoundaryConditions();
if (returnReduce(changedValues, orOp<bool>()))
{
// We've changed internal values so give
// boundary conditions opportunity to correct
he.correctBoundaryConditions();
}
}

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2017 OpenFOAM Foundation
Copyright (C) 2018-2021 OpenCFD Ltd.
Copyright (C) 2018-2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -80,41 +80,86 @@ void Foam::fv::limitVelocity::correct(volVectorField& U)
{
const scalar maxSqrU = sqr(max_);
// Count nTotCells ourselves
// (maybe only applying on a subset)
label nCellsAbove(0);
const label nTotCells(returnReduce(cells_.size(), sumOp<label>()));
vectorField& Uif = U.primitiveFieldRef();
for (const label celli : cells_)
{
const scalar magSqrUi = magSqr(Uif[celli]);
auto& Uval = Uif[celli];
const scalar magSqrUi = magSqr(Uval);
if (magSqrUi > maxSqrU)
{
Uif[celli] *= sqrt(maxSqrU/magSqrUi);
Uval *= sqrt(maxSqrU/magSqrUi);
++nCellsAbove;
}
}
// Handle boundaries in the case of 'all'
label nFacesAbove(0);
label nTotFaces(0);
if (!cellSetOption::useSubMesh())
{
for (fvPatchVectorField& Up : U.boundaryFieldRef())
{
if (!Up.fixesValue())
{
forAll(Up, facei)
// Do not count patches that fix velocity themselves
nTotFaces += Up.size();
for (auto& Uval : Up)
{
const scalar magSqrUi = magSqr(Up[facei]);
const scalar magSqrUi = magSqr(Uval);
if (magSqrUi > maxSqrU)
{
Up[facei] *= sqrt(maxSqrU/magSqrUi);
Uval *= sqrt(maxSqrU/magSqrUi);
++nFacesAbove;
}
}
}
}
}
// We've changed internal values so give
// boundary conditions opportunity to correct
U.correctBoundaryConditions();
// Percent, max 2 decimal places
const auto percent = [](scalar num, label denom) -> scalar
{
return (denom ? 1e-2*round(1e4*num/denom) : 0);
};
reduce(nCellsAbove, sumOp<label>());
// Report total numbers and percent
Info<< type() << ' ' << name_ << " Limited ";
Info<< nCellsAbove << " ("
<< percent(nCellsAbove, nTotCells)
<< "%) of cells";
reduce(nTotFaces, sumOp<label>());
reduce(nFacesAbove, sumOp<label>());
if (nTotFaces)
{
Info<< ", " << nFacesAbove << " ("
<< percent(nFacesAbove, nTotFaces)
<< "%) of faces";
}
Info<< ", with max limit " << max_ << endl;
if (nCellsAbove || nFacesAbove)
{
// We've changed internal values so give
// boundary conditions opportunity to correct
U.correctBoundaryConditions();
}
}