Skip to content
Snippets Groups Projects
Commit 248b15e7 authored by Mark Olesen's avatar Mark Olesen
Browse files

added massFlowRate bc to ease transition from OpenFOAM-1.4

parent eea37cfc
No related branches found
No related tags found
No related merge requests found
......@@ -108,6 +108,8 @@ $(derivedFvPatchFields)/turbulentIntensityKineticEnergyInlet/turbulentIntensityK
$(derivedFvPatchFields)/uniformFixedValue/uniformFixedValueFvPatchFields.C
$(derivedFvPatchFields)/waveTransmissive/waveTransmissiveFvPatchFields.C
$(derivedFvPatchFields)/massFlowRateInletVelocity/massFlowRateInletVelocityFvPatchVectorField.C
fvsPatchFields = fields/fvsPatchFields
$(fvsPatchFields)/fvsPatchField/fvsPatchFields.C
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2006-2008 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "massFlowRateInletVelocityFvPatchVectorField.H"
#include "volFields.H"
#include "addToRunTimeSelectionTable.H"
#include "fvPatchFieldMapper.H"
#include "surfaceFields.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::
massFlowRateInletVelocityFvPatchVectorField::
massFlowRateInletVelocityFvPatchVectorField
(
const fvPatch& p,
const DimensionedField<vector, volMesh>& iF
)
:
fixedValueFvPatchField<vector>(p, iF),
flowRate_(0),
phiName_("phi"),
rhoName_("rho")
{}
Foam::
massFlowRateInletVelocityFvPatchVectorField::
massFlowRateInletVelocityFvPatchVectorField
(
const massFlowRateInletVelocityFvPatchVectorField& ptf,
const fvPatch& p,
const DimensionedField<vector, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
fixedValueFvPatchField<vector>(ptf, p, iF, mapper),
flowRate_(ptf.flowRate_),
phiName_(ptf.phiName_),
rhoName_(ptf.rhoName_)
{}
Foam::
massFlowRateInletVelocityFvPatchVectorField::
massFlowRateInletVelocityFvPatchVectorField
(
const fvPatch& p,
const DimensionedField<vector, volMesh>& iF,
const dictionary& dict
)
:
fixedValueFvPatchField<vector>(p, iF, dict),
flowRate_(readScalar(dict.lookup("massFlowRate"))),
phiName_(dict.lookupOrDefault<word>("phi", "phi")),
rhoName_(dict.lookupOrDefault<word>("rho", "rho"))
{
WarningIn("massFlowRateInletVelocityFvPatchVectorField(..., dict)")
<< " changed to flowRateInletVelocityFvPatchVectorField"
<< endl;
}
Foam::
massFlowRateInletVelocityFvPatchVectorField::
massFlowRateInletVelocityFvPatchVectorField
(
const massFlowRateInletVelocityFvPatchVectorField& ptf
)
:
fixedValueFvPatchField<vector>(ptf),
flowRate_(ptf.flowRate_),
phiName_(ptf.phiName_),
rhoName_(ptf.rhoName_)
{}
Foam::
massFlowRateInletVelocityFvPatchVectorField::
massFlowRateInletVelocityFvPatchVectorField
(
const massFlowRateInletVelocityFvPatchVectorField& ptf,
const DimensionedField<vector, volMesh>& iF
)
:
fixedValueFvPatchField<vector>(ptf, iF),
flowRate_(ptf.flowRate_),
phiName_(ptf.phiName_),
rhoName_(ptf.rhoName_)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::massFlowRateInletVelocityFvPatchVectorField::updateCoeffs()
{
if (updated())
{
return;
}
// a simpler way of doing this would be nice
scalar avgU = -flowRate_/gSum(patch().magSf());
vectorField n = patch().nf();
const surfaceScalarField& phi = db().lookupObject<surfaceScalarField>
(
phiName_
);
if (phi.dimensions() == dimVelocity*dimArea)
{
// volumetric flow-rate
operator==(n*avgU);
}
else if (phi.dimensions() == dimDensity*dimVelocity*dimArea)
{
const fvPatchField<scalar>& rhop =
patch().lookupPatchField<volScalarField, scalar>(rhoName_);
// mass flow-rate
operator==(n*avgU/rhop);
}
else
{
FatalErrorIn
(
"massFlowRateInletVelocityFvPatchVectorField::updateCoeffs()"
) << "dimensions of " << phiName_ << " are incorrect" << nl
<< " on patch " << this->patch().name()
<< " of field " << this->dimensionedInternalField().name()
<< " in file " << this->dimensionedInternalField().objectPath()
<< nl << exit(FatalError);
}
fixedValueFvPatchField<vector>::updateCoeffs();
}
void Foam::massFlowRateInletVelocityFvPatchVectorField::write(Ostream& os) const
{
os.writeKeyword("type") << "flowRateInletVelocity;" << nl;
os.writeKeyword("flowRate") << flowRate_
<< token::END_STATEMENT << nl;
if (phiName_ != "phi")
{
os.writeKeyword("phi") << phiName_ << token::END_STATEMENT << nl;
}
if (rhoName_ != "rho")
{
os.writeKeyword("rho") << rhoName_ << token::END_STATEMENT << nl;
}
writeEntry("value", os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
makePatchTypeField
(
fvPatchVectorField,
massFlowRateInletVelocityFvPatchVectorField
);
}
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2006-2008 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::massFlowRateInletVelocityFvPatchVectorField
Description
Compatibility boundary condition for upgrading from
massFlowRateInletVelocity to flowRateInletVelocity
Reads 'massFlowRate' and writes 'flowRate'
Deprecated
This boundary condition is present for backward-compatibility only.
Use Foam::flowRateInletVelocityFvPatchVectorField
SourceFiles
massFlowRateInletVelocityFvPatchVectorField.C
\*---------------------------------------------------------------------------*/
#ifndef massFlowRateInletVelocityFvPatchVectorField_H
#define massFlowRateInletVelocityFvPatchVectorField_H
#include "fixedValueFvPatchFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class massFlowRateInletVelocityFvPatch Declaration
\*---------------------------------------------------------------------------*/
class massFlowRateInletVelocityFvPatchVectorField
:
public fixedValueFvPatchVectorField
{
// Private data
//- Inlet integral flow rate
scalar flowRate_;
//- Name of the flux transporting the field
word phiName_;
//- Name of the density field used to normalize the mass flux
word rhoName_;
public:
//- Runtime type information
TypeName("massFlowRateInletVelocity");
// Constructors
//- Construct from patch and internal field
massFlowRateInletVelocityFvPatchVectorField
(
const fvPatch&,
const DimensionedField<vector, volMesh>&
);
//- Construct from patch, internal field and dictionary
massFlowRateInletVelocityFvPatchVectorField
(
const fvPatch&,
const DimensionedField<vector, volMesh>&,
const dictionary&
);
//- Construct by mapping given
// massFlowRateInletVelocityFvPatchVectorField
// onto a new patch
massFlowRateInletVelocityFvPatchVectorField
(
const massFlowRateInletVelocityFvPatchVectorField&,
const fvPatch&,
const DimensionedField<vector, volMesh>&,
const fvPatchFieldMapper&
);
//- Construct as copy
massFlowRateInletVelocityFvPatchVectorField
(
const massFlowRateInletVelocityFvPatchVectorField&
);
//- Construct and return a clone
virtual tmp<fvPatchVectorField> clone() const
{
return tmp<fvPatchVectorField>
(
new massFlowRateInletVelocityFvPatchVectorField(*this)
);
}
//- Construct as copy setting internal field reference
massFlowRateInletVelocityFvPatchVectorField
(
const massFlowRateInletVelocityFvPatchVectorField&,
const DimensionedField<vector, volMesh>&
);
//- Construct and return a clone setting internal field reference
virtual tmp<fvPatchVectorField> clone
(
const DimensionedField<vector, volMesh>& iF
) const
{
return tmp<fvPatchVectorField>
(
new massFlowRateInletVelocityFvPatchVectorField(*this, iF)
);
}
// Member functions
// Access
//- Return the flux
scalar flowRate() const
{
return flowRate_;
}
//- Return reference to the flux to allow adjustment
scalar& flowRate()
{
return flowRate_;
}
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();
//- Write
virtual void write(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment