/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see .
\*---------------------------------------------------------------------------*/
#include "dsmcFields.H"
#include "volFields.H"
#include "dictionary.H"
#include "dsmcCloud.H"
#include "constants.H"
#include "addToRunTimeSelectionTable.H"
using namespace Foam::constant;
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace functionObjects
{
defineTypeNameAndDebug(dsmcFields, 0);
addToRunTimeSelectionTable
(
functionObject,
dsmcFields,
dictionary
);
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::functionObjects::dsmcFields::dsmcFields
(
const word& name,
const Time& runTime,
const dictionary& dict
)
:
functionObject(name),
obr_
(
runTime.lookupObject
(
dict.lookupOrDefault("region", polyMesh::defaultRegion)
)
)
{
if (!isA(obr_))
{
FatalErrorInFunction
<< "objectRegistry is not an fvMesh" << exit(FatalError);
}
read(dict);
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::functionObjects::dsmcFields::~dsmcFields()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::functionObjects::dsmcFields::read(const dictionary& dict)
{
return true;
}
bool Foam::functionObjects::dsmcFields::execute()
{
return true;
}
bool Foam::functionObjects::dsmcFields::write()
{
word rhoNName = "rhoN";
word rhoMName = "rhoM";
word momentumName = "momentum";
word linearKEName = "linearKE";
word internalEName = "internalE";
word iDofName = "iDof";
word fDName = "fD";
//obr_.lookupObject()
/*volScalarField rhoNMean
(
IOobject
(
rhoNMeanName,
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1
);*/
const volScalarField& baseField1 = obr_.lookupObject(rhoNName);
obr_.store
(
new volScalarField
(
IOobject
(
rhoNName + "Mean",
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1*baseField1
)
);
const volScalarField& baseField2 = obr_.lookupObject(rhoMName);
obr_.store
(
new volScalarField
(
IOobject
(
rhoMName + "Mean",
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1*baseField2
)
);
const volVectorField& baseField3 = obr_.lookupObject(momentumName);
obr_.store
(
new volVectorField
(
IOobject
(
momentumName + "Mean",
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1*baseField3
)
);
const volScalarField& baseField4 = obr_.lookupObject(linearKEName);
obr_.store
(
new volScalarField
(
IOobject
(
linearKEName + "Mean",
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1*baseField4
)
);
const volScalarField& baseField5 = obr_.lookupObject(internalEName);
obr_.store
(
new volScalarField
(
IOobject
(
internalEName + "Mean",
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1*baseField5
)
);
const volScalarField& baseField6 = obr_.lookupObject(iDofName);
obr_.store
(
new volScalarField
(
IOobject
(
iDofName + "Mean",
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1*baseField6
)
);
const volVectorField& baseField7 = obr_.lookupObject(fDName);
obr_.store
(
new volVectorField
(
IOobject
(
fDName + "Mean",
obr_.time().timeName(obr_.time().startTime().value()),
obr_,
IOobject::MUST_READ,
IOobject::NO_WRITE
),
1*baseField7
)
);
const volScalarField& rhoNMean = obr_.lookupObject
(
rhoNName + "Mean"
);
const volScalarField& rhoMMean = obr_.lookupObject
(
rhoMName + "Mean"
);
const volVectorField& momentumMean = obr_.lookupObject
(
momentumName + "Mean"
);
const volScalarField& linearKEMean = obr_.lookupObject
(
linearKEName + "Mean"
);
const volScalarField& internalEMean = obr_.lookupObject
(
internalEName + "Mean"
);
const volScalarField& iDofMean = obr_.lookupObject
(
iDofName + "Mean"
);
const volVectorField& fDMean = obr_.lookupObject
(
fDName + "Mean"
);
if (min(mag(rhoNMean)).value() > VSMALL)
{
Info<< "Calculating dsmcFields." << endl;
Info<< " Calculating UMean field." << endl;
volVectorField UMean
(
IOobject
(
"UMean",
obr_.time().timeName(),
obr_,
IOobject::NO_READ
),
momentumMean/rhoMMean
);
Info<< " Calculating translationalT field." << endl;
volScalarField translationalT
(
IOobject
(
"translationalT",
obr_.time().timeName(),
obr_,
IOobject::NO_READ
),
2.0/(3.0*physicoChemical::k.value()*rhoNMean)
*(linearKEMean - 0.5*rhoMMean*(UMean & UMean))
);
Info<< " Calculating internalT field." << endl;
volScalarField internalT
(
IOobject
(
"internalT",
obr_.time().timeName(),
obr_,
IOobject::NO_READ
),
(2.0/physicoChemical::k.value())*(internalEMean/iDofMean)
);
Info<< " Calculating overallT field." << endl;
volScalarField overallT
(
IOobject
(
"overallT",
obr_.time().timeName(),
obr_,
IOobject::NO_READ
),
2.0/(physicoChemical::k.value()*(3.0*rhoNMean + iDofMean))
*(linearKEMean - 0.5*rhoMMean*(UMean & UMean) + internalEMean)
);
Info<< " Calculating pressure field." << endl;
volScalarField p
(
IOobject
(
"p",
obr_.time().timeName(),
obr_,
IOobject::NO_READ
),
physicoChemical::k.value()*rhoNMean*translationalT
);
const fvMesh& mesh = fDMean.mesh();
volScalarField::Boundary& pBf = p.boundaryFieldRef();
forAll(mesh.boundaryMesh(), i)
{
const polyPatch& patch = mesh.boundaryMesh()[i];
if (isA(patch))
{
pBf[i] =
fDMean.boundaryField()[i]
& (patch.faceAreas()/mag(patch.faceAreas()));
}
}
Info<< " mag(UMean) max/min : "
<< max(mag(UMean)).value() << " "
<< min(mag(UMean)).value() << endl;
Info<< " translationalT max/min : "
<< max(translationalT).value() << " "
<< min(translationalT).value() << endl;
Info<< " internalT max/min : "
<< max(internalT).value() << " "
<< min(internalT).value() << endl;
Info<< " overallT max/min : "
<< max(overallT).value() << " "
<< min(overallT).value() << endl;
Info<< " p max/min : "
<< max(p).value() << " "
<< min(p).value() << endl;
UMean.write();
translationalT.write();
internalT.write();
overallT.write();
p.write();
Info<< "dsmcFields written." << nl << endl;
return true;
}
else
{
Info<< "Small value (" << min(mag(rhoNMean))
<< ") found in rhoNMean field. "
<< "Not calculating dsmcFields to avoid division by zero."
<< endl;
return false;
}
}
// ************************************************************************* //