Skip to content
Permalink
Browse files

functionObjects: Added age function object

This object calculates a field of the age of fluid in the domain; i.e.,
the time taken for a fluid particle to travel to a location from an
inlet. It outputs a field, named age, with dimensions of time, and
requires a solver and a div(phi,age) scheme to be specified. A number of
corrections for the solution procedure can be set, as well as the name
of the flux and density fields.

Example specification:

    age1
    {
        type    age;
        libs    ("libfieldFunctionObjects.so");
        nCorr   10;
        phi     phi;
        rho     rho;
    }

Example usage:

    postProcess -func age -fields "(phi)" -latestTime

This work was supported by Robert Secor and Lori Holmes, at 3M
  • Loading branch information...
Will Bainbridge
Will Bainbridge committed Nov 20, 2018
1 parent e928eaf commit b928e3767725ece80c48c4970715d94a51907231
@@ -0,0 +1,26 @@
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: dev
\\/ M anipulation |
-------------------------------------------------------------------------------
Description
Solves a transport equation to determine the time taken for a particle to
convect from an inlet to the location in the flow.
This will output a field, age, with units of time. This field needs a
solver setting in fvSolution and a div(phi,age) scheme in fvSchemes.
The number of correctors, nCorr, determines how many times the solution is
repeated to iterate away any non-linearities in the choice of scheme. If
the divergence scheme is set to upwind, no corrections will be necessary.
\*---------------------------------------------------------------------------*/

type age;
libs ("libfieldFunctionObjects.so");

nCorr 10;

// ************************************************************************* //
@@ -67,4 +67,6 @@ subtract/subtract.C

interfaceHeight/interfaceHeight.C

age/age.C

LIB = $(FOAM_LIBBIN)/libfieldFunctionObjects
@@ -0,0 +1,162 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2018 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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/

#include "age.H"
#include "fvCFD.H"
#include "inletOutletFvPatchField.H"
#include "wallFvPatch.H"
#include "zeroGradientFvPatchField.H"
#include "addToRunTimeSelectionTable.H"

// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

namespace Foam
{
namespace functionObjects
{
defineTypeNameAndDebug(age, 0);
addToRunTimeSelectionTable(functionObject, age, dictionary);
}
}


// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //

Foam::wordList Foam::functionObjects::age::patchTypes() const
{
wordList result
(
mesh_.boundary().size(),
inletOutletFvPatchField<scalar>::typeName
);

forAll(mesh_.boundary(), patchi)
{
if (isA<wallFvPatch>(mesh_.boundary()[patchi]))
{
result[patchi] = zeroGradientFvPatchField<scalar>::typeName;
}
}

return result;
}


// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //

Foam::functionObjects::age::age
(
const word& name,
const Time& runTime,
const dictionary& dict
)
:
fvMeshFunctionObject(name, runTime, dict),
nCorr_(readLabel(dict.lookup("nCorr"))),
phiName_(),
rhoName_()
{
read(dict);
}


// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //

Foam::functionObjects::age::~age()
{}


// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //

bool Foam::functionObjects::age::read(const dictionary& dict)
{
phiName_ = dict.lookupOrDefault<word>("phi", "phi");
rhoName_ = dict.lookupOrDefault<word>("rho", "rho");

dict.readIfPresent("nCorr", nCorr_);

return true;
}


bool Foam::functionObjects::age::execute()
{
volScalarField t
(
IOobject
(
"age",
mesh_.time().timeName(),
mesh_,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh_,
dimensionedScalar("zero", dimTime, 0),
patchTypes()
);

// This only works because the null constructed inletValue for an
// inletOutletFvPatchField is zero. If we needed any other value we would
// have to loop over the inletOutlet patches and explicitly set the
// inletValues. We would need to change the interface of inletOutlet in
// order to do this.

const surfaceScalarField& phi =
mesh_.lookupObject<surfaceScalarField>(phiName_);

if (phi.dimensions() == dimMass/dimTime)
{
const volScalarField& rho =
mesh_.lookupObject<volScalarField>(rhoName_);

for (label i = 0; i <= nCorr_; ++ i)
{
solve(fvm::div(phi, t) == rho);
}
}
else
{
for (label i = 0; i <= nCorr_; ++ i)
{
solve(fvm::div(phi, t) == dimensionedScalar("one", dimless, 1));
}
}

Info<< "Min/max age:" << min(t).value() << ' ' << max(t).value() << endl;

t.write();

return true;
}


bool Foam::functionObjects::age::write()
{
return true;
}


// ************************************************************************* //
@@ -0,0 +1,151 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2018 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 <http://www.gnu.org/licenses/>.
Class
Foam::functionObjects::age
Description
Calculates and writes out the time taken for a particle to travel from an
inlet to the location. Solves the following equation when incompressible:
\f[
\div (\phi t) = 1
\f]
where:
\vartable
t | Age [s]
\phi | Volumetric flux [m^3/s]
\endvartable
Boundary conditions are generated automatically as zeroGradient on all
walls and inletOutlet everywhere else.
Usage
\table
Property | Description | Required | Default value
nCorr | The number of correctors | yes |
phi | The name of the flux field | no | phi
rho | The name of the density field | no | rho
\endtable
\verbatim
age1
{
type age;
libs ("libsolverFunctionObjects.so");
writeControl writeTime;
writeInterval 1;
nCorr 10;
phi phi;
rho rho;
}
\endverbatim
SourceFiles
age.C
\*---------------------------------------------------------------------------*/

#ifndef age_H
#define age_H

#include "fvMeshFunctionObject.H"
#include "volFields.H"

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

namespace Foam
{
namespace functionObjects
{

/*---------------------------------------------------------------------------*\
Class age Declaration
\*---------------------------------------------------------------------------*/

class age
:
public fvMeshFunctionObject
{
// Private data

//- Number of corrections
label nCorr_;

//- The name of the flux field
word phiName_;

//- The name of the density field
word rhoName_;


// Private Member Functions

//- The list of patch types for the age field
wordList patchTypes() const;


public:

//- Runtime type information
TypeName("age");


// Constructors

//- Construct from Time and dictionary
age
(
const word& name,
const Time& runTime,
const dictionary& dict
);


//- Destructor
virtual ~age();


// Member Functions

//- Read the data
virtual bool read(const dictionary&);

//- Execute
virtual bool execute();

//- Write
virtual bool write();
};


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

} // End namespace functionObjects
} // End namespace Foam

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

#endif

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

0 comments on commit b928e37

Please sign in to comment.
You can’t perform that action at this time.