Commit 663e9b7a authored by Mark Olesen's avatar Mark Olesen
Browse files

start simplification of fieldAverage

- track previous time index and to avoid recalculating averages.
- make sure averages are up-to-date before calling write
parent fee6e312
......@@ -26,19 +26,13 @@ License
#include "fieldAverage.H"
#include "volFields.H"
#include "dictionary.H"
#include "Time.H"
#include "IFstream.H"
#include "OFstream.H"
#include "fieldAverageItem.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(fieldAverage, 0);
}
defineTypeNameAndDebug(Foam::fieldAverage, 0);
const Foam::word Foam::fieldAverage::EXT_MEAN = "Mean";
const Foam::word Foam::fieldAverage::EXT_PRIME2MEAN = "Prime2Mean";
......@@ -46,100 +40,81 @@ const Foam::word Foam::fieldAverage::EXT_PRIME2MEAN = "Prime2Mean";
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::fieldAverage::checkoutFields(const wordList& fieldNames) const
void Foam::fieldAverage::resetFields(wordList& names)
{
forAll(fieldNames, i)
forAll(names, fieldI)
{
if (fieldNames[i] != word::null)
if (names[fieldI].size())
{
obr_.checkOut(*obr_[fieldNames[i]]);
obr_.checkOut(*obr_[names[fieldI]]);
}
}
names.clear();
names.setSize(faItems_.size());
}
void Foam::fieldAverage::resetLists(const label nItems)
void Foam::fieldAverage::initialize()
{
checkoutFields(meanScalarFields_);
meanScalarFields_.clear();
meanScalarFields_.setSize(nItems);
checkoutFields(meanVectorFields_);
meanVectorFields_.clear();
meanVectorFields_.setSize(nItems);
checkoutFields(meanSphericalTensorFields_);
meanSphericalTensorFields_.clear();
meanSphericalTensorFields_.setSize(nItems);
checkoutFields(meanSymmTensorFields_);
meanSymmTensorFields_.clear();
meanSymmTensorFields_.setSize(nItems);
resetFields(meanScalarFields_);
resetFields(meanVectorFields_);
resetFields(meanSphericalTensorFields_);
resetFields(meanSymmTensorFields_);
resetFields(meanTensorFields_);
checkoutFields(meanTensorFields_);
meanTensorFields_.clear();
meanTensorFields_.setSize(nItems);
checkoutFields(prime2MeanScalarFields_);
prime2MeanScalarFields_.clear();
prime2MeanScalarFields_.setSize(nItems);
checkoutFields(prime2MeanSymmTensorFields_);
prime2MeanSymmTensorFields_.clear();
prime2MeanSymmTensorFields_.setSize(nItems);
resetFields(prime2MeanScalarFields_);
resetFields(prime2MeanSymmTensorFields_);
totalIter_.clear();
totalIter_.setSize(nItems, 1);
totalIter_.setSize(faItems_.size(), 1);
totalTime_.clear();
totalTime_.setSize(nItems, obr_.time().deltaT().value());
}
totalTime_.setSize(faItems_.size(), obr_.time().deltaT().value());
void Foam::fieldAverage::initialise()
{
// Add mean fields to the field lists
forAll(faItems_, i)
forAll(faItems_, fieldI)
{
const word& fieldName = faItems_[i].fieldName();
const word& fieldName = faItems_[fieldI].fieldName();
if (obr_.foundObject<volScalarField>(fieldName))
{
addMeanField<scalar>(i, meanScalarFields_);
addMeanField<scalar>(fieldI, meanScalarFields_);
}
else if (obr_.foundObject<volVectorField>(fieldName))
{
addMeanField<vector>(i, meanVectorFields_);
addMeanField<vector>(fieldI, meanVectorFields_);
}
else if (obr_.foundObject<volSphericalTensorField>(fieldName))
{
addMeanField<sphericalTensor>(i, meanSphericalTensorFields_);
addMeanField<sphericalTensor>(fieldI, meanSphericalTensorFields_);
}
else if (obr_.foundObject<volSymmTensorField>(fieldName))
{
addMeanField<symmTensor>(i, meanSymmTensorFields_);
addMeanField<symmTensor>(fieldI, meanSymmTensorFields_);
}
else if (obr_.foundObject<volTensorField>(fieldName))
{
addMeanField<tensor>(i, meanTensorFields_);
addMeanField<tensor>(fieldI, meanTensorFields_);
}
else
{
FatalErrorIn("Foam::fieldAverage::initialise()")
<< "Requested field " << faItems_[i].fieldName()
FatalErrorIn("Foam::fieldAverage::initialize()")
<< "Requested field " << faItems_[fieldI].fieldName()
<< " does not exist in the database" << nl
<< exit(FatalError);
}
}
// Add prime-squared mean fields to the field lists
forAll(faItems_, i)
forAll(faItems_, fieldI)
{
if (faItems_[i].prime2Mean())
if (faItems_[fieldI].prime2Mean())
{
const word& fieldName = faItems_[i].fieldName();
if (!faItems_[i].mean())
const word& fieldName = faItems_[fieldI].fieldName();
if (!faItems_[fieldI].mean())
{
FatalErrorIn("Foam::fieldAverage::initialise()")
FatalErrorIn("Foam::fieldAverage::initialize()")
<< "To calculate the prime-squared average, the "
<< "mean average must also be selected for field "
<< fieldName << nl << exit(FatalError);
......@@ -149,7 +124,7 @@ void Foam::fieldAverage::initialise()
{
addPrime2MeanField<scalar, scalar>
(
i,
fieldI,
meanScalarFields_,
prime2MeanScalarFields_
);
......@@ -158,14 +133,14 @@ void Foam::fieldAverage::initialise()
{
addPrime2MeanField<vector, symmTensor>
(
i,
fieldI,
meanVectorFields_,
prime2MeanSymmTensorFields_
);
}
else
{
FatalErrorIn("Foam::fieldAverage::initialise()")
FatalErrorIn("Foam::fieldAverage::initialize()")
<< "prime2Mean average can only be applied to "
<< "volScalarFields and volVectorFields"
<< nl << " Field: " << fieldName << nl
......@@ -176,111 +151,26 @@ void Foam::fieldAverage::initialise()
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::fieldAverage::fieldAverage
(
const word& name,
const objectRegistry& obr,
const dictionary& dict,
const bool loadFromFiles
)
:
name_(name),
obr_(obr),
active_(true),
cleanRestart_(dict.lookupOrDefault<Switch>("cleanRestart", false)),
faItems_(dict.lookup("fields")),
meanScalarFields_(faItems_.size()),
meanVectorFields_(faItems_.size()),
meanSphericalTensorFields_(faItems_.size()),
meanSymmTensorFields_(faItems_.size()),
meanTensorFields_(faItems_.size()),
prime2MeanScalarFields_(faItems_.size()),
prime2MeanSymmTensorFields_(faItems_.size()),
totalIter_(faItems_.size(), 1),
totalTime_(faItems_.size(), obr_.time().deltaT().value())
{
// Check if the available mesh is an fvMesh otherise deactivate
if (!isA<fvMesh>(obr_))
{
active_ = false;
WarningIn
(
"fieldAverage::fieldAverage\n"
"(\n"
"const word&,\n"
"const objectRegistry&,\n"
"const dictionary&,\n"
"const bool\n"
")"
) << "No fvMesh available, deactivating."
<< nl << endl;
}
read(dict);
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::fieldAverage::~fieldAverage()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::fieldAverage::read(const dictionary& dict)
{
if (active_)
{
faItems_.clear();
faItems_ = List<fieldAverageItem>(dict.lookup("fields"));
resetLists(faItems_.size());
initialise();
readAveragingProperties();
}
}
void Foam::fieldAverage::execute()
void Foam::fieldAverage::calcAverages()
{
if (active_)
{
calcAverages();
}
}
const label currentTimeIndex =
static_cast<const fvMesh&>(obr_).time().timeIndex();
void Foam::fieldAverage::end()
{
if (active_)
if (prevTimeIndex_ == currentTimeIndex)
{
calcAverages();
return;
}
}
void Foam::fieldAverage::write()
{
if (active_)
else
{
writeAverages();
writeAveragingProperties();
prevTimeIndex_ = currentTimeIndex;
}
}
void Foam::fieldAverage::calcAverages()
{
Info<< "Calculating averages" << nl << endl;
forAll(faItems_, i)
forAll(faItems_, fieldI)
{
totalIter_[i]++;
totalTime_[i] += obr_.time().deltaT().value();
totalIter_[fieldI]++;
totalTime_[fieldI] += obr_.time().deltaT().value();
}
addMeanSqrToPrime2Mean<scalar, scalar>
......@@ -342,12 +232,12 @@ void Foam::fieldAverage::writeAveragingProperties() const
)
);
forAll(faItems_, i)
forAll(faItems_, fieldI)
{
const word& fieldName = faItems_[i].fieldName();
const word& fieldName = faItems_[fieldI].fieldName();
propsDict.add(fieldName, dictionary());
propsDict.subDict(fieldName).add("totalIter", totalIter_[i]);
propsDict.subDict(fieldName).add("totalTime", totalTime_[i]);
propsDict.subDict(fieldName).add("totalIter", totalIter_[fieldI]);
propsDict.subDict(fieldName).add("totalTime", totalTime_[fieldI]);
}
propsDict.regIOobject::write();
......@@ -363,34 +253,39 @@ void Foam::fieldAverage::readAveragingProperties()
}
else
{
IFstream propsFile
IOobject propsDictHeader
(
obr_.time().path()/obr_.time().timeName()
/"uniform"/"fieldAveragingProperties"
"fieldAveragingProperties",
obr_.time().timeName(),
"uniform",
obr_,
IOobject::NO_READ,
IOobject::NO_WRITE,
false
);
if (!propsFile.good())
if (!propsDictHeader.headerOk())
{
Info<< "fieldAverage: starting averaging at time "
<< obr_.time().timeName() << nl << endl;
return;
}
dictionary propsDict(dictionary::null, propsFile);
IOdictionary propsDict(propsDictHeader);
Info<< "fieldAverage: restarting averaging for fields:" << endl;
forAll(faItems_, i)
forAll(faItems_, fieldI)
{
const word& fieldName = faItems_[i].fieldName();
const word& fieldName = faItems_[fieldI].fieldName();
if (propsDict.found(fieldName))
{
dictionary fieldDict(propsDict.subDict(fieldName));
totalIter_[i] = readLabel(fieldDict.lookup("totalIter"));
totalTime_[i] = readScalar(fieldDict.lookup("totalTime"));
totalIter_[fieldI] = readLabel(fieldDict.lookup("totalIter"));
totalTime_[fieldI] = readScalar(fieldDict.lookup("totalTime"));
Info<< " " << fieldName
<< " iters = " << totalIter_[i]
<< " time = " << totalTime_[i] << endl;
<< " iters = " << totalIter_[fieldI]
<< " time = " << totalTime_[fieldI] << endl;
}
}
Info<< endl;
......@@ -398,6 +293,104 @@ void Foam::fieldAverage::readAveragingProperties()
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::fieldAverage::fieldAverage
(
const word& name,
const objectRegistry& obr,
const dictionary& dict,
const bool loadFromFiles
)
:
name_(name),
obr_(obr),
active_(true),
prevTimeIndex_(-1),
cleanRestart_(false),
faItems_(),
meanScalarFields_(),
meanVectorFields_(),
meanSphericalTensorFields_(),
meanSymmTensorFields_(),
meanTensorFields_(),
prime2MeanScalarFields_(),
prime2MeanSymmTensorFields_(),
totalIter_(),
totalTime_()
{
// Only active if a fvMesh is available
if (isA<fvMesh>(obr_))
{
read(dict);
}
else
{
active_ = false;
WarningIn
(
"fieldAverage::fieldAverage\n"
"(\n"
"const word&,\n"
"const objectRegistry&,\n"
"const dictionary&,\n"
"const bool\n"
")"
) << "No fvMesh available, deactivating."
<< nl << endl;
}
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::fieldAverage::~fieldAverage()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::fieldAverage::read(const dictionary& dict)
{
if (active_)
{
dict.readIfPresent("cleanRestart", cleanRestart_);
dict.lookup("fields") >> faItems_;
initialize();
readAveragingProperties();
// ensure first averaging works unconditionally
prevTimeIndex_ = -1;
}
}
void Foam::fieldAverage::execute()
{
if (active_)
{
calcAverages();
}
}
void Foam::fieldAverage::end()
{
}
void Foam::fieldAverage::write()
{
if (active_)
{
calcAverages();
writeAverages();
writeAveragingProperties();
}
}
void Foam::fieldAverage::updateMesh(const mapPolyMesh&)
{
// Do nothing
......
......@@ -107,6 +107,14 @@ class fieldAverage
{
protected:
// File and field name extensions
//- Mean average
static const word EXT_MEAN;
//- Prime-squared average
static const word EXT_PRIME2MEAN;
// Private data
//- Name of this set of field averages.
......@@ -118,6 +126,9 @@ protected:
//- On/off switch
bool active_;
//- Time at last call, prevents repeated averaging
label prevTimeIndex_;
//- Clean restart flag
Switch cleanRestart_;
......@@ -125,15 +136,6 @@ protected:
// calculated and output
List<fieldAverageItem> faItems_;
// File and field name extensions
//- Mean average
static const word EXT_MEAN;
//- Prime-squared average
static const word EXT_PRIME2MEAN;
// Lists of averages
// Arithmetic mean fields
......@@ -143,7 +145,8 @@ protected:
wordList meanSymmTensorFields_;
wordList meanTensorFields_;
// Prime-squared fields - applicable to volVectorFields only
// Prime-squared fields
// Only applicable to volScalarFields / volVectorFields
wordList prime2MeanScalarFields_;
wordList prime2MeanSymmTensorFields_;
......@@ -162,20 +165,18 @@ protected:
// Initialisation routines
//- Checkout fields (causes deletion) from the database
void checkoutFields(const wordList&) const;
//- Reset size of lists (clear existing values)
void resetLists(const label nItems);
// and reset lists
void resetFields(wordList&);
//- Intitialise averaging. Check requested field averages are
// valid, and populate field lists
void initialise();
//- Reset lists (clear existing values) and initialize averaging.
// Check requested field averages are valid, populate field lists
void initialize();
//- Add mean average field to PtrList
//- Add mean average field to list
template<class Type>
void addMeanField(const label, wordList&) const;
//- Add prime-squared average field to PtrList
//- Add prime-squared average field to list
template<class Type1, class Type2>
void addPrime2MeanField
(
......@@ -211,7 +212,7 @@ protected:
) const;
// I-O
// IO
//- Write averages
virtual void writeAverages() const;
......@@ -281,7 +282,7 @@ public:
//- Execute the averaging
virtual void execute();
//- Execute the averaging at the final time-loop
//- Execute the averaging at the final time-loop, currently does nothing
virtual void end();
//- Calculate the field average data and write
......
......@@ -33,15 +33,15 @@ License
template<class Type>
void Foam::fieldAverage::addMeanField
(
const label fieldi,
const label fieldI,
wordList& meanFieldList
) const
{
if (faItems_[fieldi].mean())
if (faItems_[fieldI].mean())
{
typedef GeometricField<Type, fvPatchField, volMesh> fieldType;
const word& fieldName = faItems_[fieldi].fieldName();
const word& fieldName = faItems_[fieldI].fieldName();
const word meanFieldName = fieldName + EXT_MEAN;
......@@ -49,14 +49,14 @@ void Foam::fieldAverage::addMeanField
if (obr_.foundObject<fieldType>(meanFieldName))
{