zhyang-dev/OpenFOAM-6
1
0
Fork 0
mirror of https://github.com/OpenFOAM/OpenFOAM-6.git synced 2025-12-08 06:57:46 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

waves: Control additions

Browse Source 
An "inletOutlet" switch has been added to the wave velocity boundary
condition to allow the boundary to be fixed, as is possible for the
corresponding alpha condition.

A "heightAboveWave" option has been added to the wave superposition
class to calculate velocity based on the height above the wave, rather
than above the origin. This may improve initialisation but it may also
generate divergence in the initial velocity field.

The alpha condition has also been completed so that it applies a
modelled gradient when the flow points out and a wave pressure condition
is in use.
This commit is contained in:
Will Bainbridge
2018-04-19 10:45:12 +01:00
parent a0a19c7f38
commit 8144cfb516
8 changed files with 217 additions and 74 deletions
Download Patch File Download Diff File Expand all files Collapse all files

91
src/waves/derivedFvPatchFields/waveAlpha/waveAlphaFvPatchScalarField.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -24,11 +24,13 @@ License
\*---------------------------------------------------------------------------*/
#include "waveAlphaFvPatchScalarField.H"
#include "wavePressureFvPatchScalarField.H"
#include "waveVelocityFvPatchVectorField.H"
#include "addToRunTimeSelectionTable.H"
#include "levelSet.H"
#include "surfaceFields.H"
#include "volFields.H"
#include "fvMeshSubset.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
@ -139,6 +141,50 @@ Foam::tmp<Foam::scalarField> Foam::waveAlphaFvPatchScalarField::alpha() const
}
Foam::tmp<Foam::scalarField> Foam::waveAlphaFvPatchScalarField::alphan() const
{
const waveVelocityFvPatchVectorField& Up =
refCast<const waveVelocityFvPatchVectorField>
(
patch().lookupPatchField<volVectorField, scalar>(UName_)
);
const scalar t = db().time().timeOutputValue();
const fvMeshSubset& subset = Up.faceCellSubset();
const fvMesh& meshs = subset.subMesh();
const label patchis = findIndex(subset.patchMap(), patch().index());
const scalarField alphas
(
levelSetFraction
(
meshs,
Up.waves().height(t, meshs.cellCentres())(),
Up.waves().height(t, meshs.points())(),
!liquid_
)
);
tmp<scalarField> tResult(new scalarField(patch().size()));
scalarField& result = tResult.ref();
if (patchis != -1)
{
forAll(meshs.boundary()[patchis], is)
{
const label fs = is + meshs.boundary()[patchis].patch().start();
const label cs = meshs.boundary()[patchis].faceCells()[is];
const label f = subset.faceMap()[fs];
const label i = patch().patch().whichFace(f);
result[i] = alphas[cs];
}
}
return tResult;
}
void Foam::waveAlphaFvPatchScalarField::updateCoeffs()
{
if (updated())
@ -146,18 +192,49 @@ void Foam::waveAlphaFvPatchScalarField::updateCoeffs()
return;
}
refValue() = alpha();
const fvPatchVectorField& Up =
patch().lookupPatchField<volVectorField, scalar>(UName_);
if (inletOutlet_)
if (!isA<waveVelocityFvPatchVectorField>(Up))
{
const scalarField& phip =
patch().lookupPatchField<surfaceScalarField, scalar>("phi");
FatalErrorInFunction
<< "The corresponding condition for the velocity "
<< "field " << UName_ << " on patch " << patch().name()
<< " is not of type " << waveVelocityFvPatchVectorField::typeName
<< exit(FatalError);
}
valueFraction() = 1 - pos0(phip);
const waveVelocityFvPatchVectorField& Uwp =
refCast<const waveVelocityFvPatchVectorField>(Up);
const fvPatchScalarField& pp =
patch().lookupPatchField<volScalarField, scalar>(Uwp.pName());
if (isA<wavePressureFvPatchScalarField>(pp))
{
const scalarField alpha(this->alpha()), alphan(this->alphan());
const scalarField out(pos0(Uwp.U() & patch().Sf()));
valueFraction() = out;
refValue() = alpha;
refGrad() = (alpha - alphan)*patch().deltaCoeffs();
}
else
{
valueFraction() = 1;
refValue() = alpha();
if (inletOutlet_)
{
const scalarField& phip =
patch().lookupPatchField<surfaceScalarField, scalar>("phi");
const scalarField out(pos0(phip));
valueFraction() = 1 - out;
}
else
{
valueFraction() = 1;
}
}
mixedFvPatchScalarField::updateCoeffs();

29
src/waves/derivedFvPatchFields/waveAlpha/waveAlphaFvPatchScalarField.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -32,6 +32,21 @@ Description
fraction to that specified by a superposition of wave models. All the
parameters are looked up from the corresponding velocity condition.
Flow reversal will occur in the event that the amplitude of the velocity
oscillation is greater than the mean flow. This triggers special handling,
the form of which depends on the inletOutlet flag and whether a wave
pressure condition is being used.
If a wave pressure condition is not being used, the inletOutlet switches
between a fixedValue and an inletOutlet condition, with the value given by
the wave model. If fixedValue, the result may be more accurate, but it
might also be unstable.
If a wave pressure condition is being used, then the normal phase fraction
condition becomes fixedGradient on outlet faces. This gradient is
calculated numerically by evaluating the wave model on both the patch face
and the adjacent cell.
Usage
\table
Property | Description | Req'd? | Default
@ -76,16 +91,12 @@ class waveAlphaFvPatchScalarField
// Private data
//- Name of the velocity field
word UName_;
const word UName_;
//- Is this alpha field that of the liquid under the wave?
const Switch liquid_;
//- Act as an inlet/outlet patch? If false, the alpha field will be set
// by the wave model even during outflow. This may apply the wave model
// more accurately, but it might also be unstable. If true, the alpha
// boundary condition will be zero gradient when the flow reverses, as
// is usual.
//- Act as an inlet/outlet patch?
const Switch inletOutlet_;
@ -172,6 +183,10 @@ public:
//- Return the current modelled phase fraction field
tmp<scalarField> alpha() const;
//- Return the current modelled phase fraction field in the
// neighbour cell
tmp<scalarField> alphan() const;
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();

3
src/waves/derivedFvPatchFields/wavePressure/wavePressureFvPatchScalarField.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -216,7 +216,6 @@ void Foam::wavePressureFvPatchScalarField::updateCoeffs()
}
const scalarField p(this->p()), pn(this->pn());
const scalarField out(pos0(Uwp.U() & patch().Sf()));
valueFraction() = out;

11
src/waves/derivedFvPatchFields/wavePressure/wavePressureFvPatchScalarField.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -32,6 +32,15 @@ Description
pressure to a value specified by a superposition of wave models. All the
parameters are looked up from the corresponding velocity condition.
This functions like an outletInlet condition. Faces on which the flow is
leaving the domain have a value set by the wave model. Faces on which the
flow is entering the domain have the gradient set. This gradient is
calculated numerically by evaluating the wave model on both the patch face
and the adjacent cell.
Use of this boundary condition triggers a consistent behaviour in the
corresponding velocity and phase-fraction conditions.
Usage
\table
Property | Description | Req'd? | Default

77
src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.C
View File

@ -40,7 +40,8 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField
:
directionMixedFvPatchVectorField(p, iF),
phiName_("phi"),
pName_(word::null),
pName_("p"),
inletOutlet_(true),
waves_(db()),
faceCellSubset_(nullptr),
faceCellSubsetTimeIndex_(-1)
@ -60,7 +61,8 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField
:
directionMixedFvPatchVectorField(p, iF),
phiName_(dict.lookupOrDefault<word>("phi", "phi")),
pName_(dict.lookupOrDefault<word>("p", word::null)),
pName_(dict.lookupOrDefault<word>("p", "p")),
inletOutlet_(dict.lookupOrDefault<Switch>("inletOutlet", true)),
waves_(db(), dict),
faceCellSubset_(nullptr),
faceCellSubsetTimeIndex_(-1)
@ -91,6 +93,7 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField
directionMixedFvPatchVectorField(ptf, p, iF, mapper),
phiName_(ptf.phiName_),
pName_(ptf.pName_),
inletOutlet_(ptf.inletOutlet_),
waves_(ptf.waves_),
faceCellSubset_(nullptr),
faceCellSubsetTimeIndex_(-1)
@ -105,6 +108,7 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField
directionMixedFvPatchVectorField(ptf),
phiName_(ptf.phiName_),
pName_(ptf.pName_),
inletOutlet_(ptf.inletOutlet_),
waves_(ptf.waves_),
faceCellSubset_(nullptr),
faceCellSubsetTimeIndex_(-1)
@ -120,6 +124,7 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField
directionMixedFvPatchVectorField(ptf, iF),
phiName_(ptf.phiName_),
pName_(ptf.pName_),
inletOutlet_(ptf.inletOutlet_),
waves_(ptf.waves_),
faceCellSubset_(nullptr),
faceCellSubsetTimeIndex_(-1)
@ -220,21 +225,11 @@ void Foam::waveVelocityFvPatchVectorField::updateCoeffs()
return;
}
if (pName_ != word::null)
const fvPatchScalarField& pp =
patch().lookupPatchField<volScalarField, scalar>(pName_);
if (isA<wavePressureFvPatchScalarField>(pp))
{
const fvPatchScalarField& pp =
patch().lookupPatchField<volScalarField, scalar>(pName_);
if (!isA<wavePressureFvPatchScalarField>(pp))
{
FatalErrorInFunction
<< "The corresponding pressure condition for the pressure "
<< "field " << pName_ << " on patch " << patch().name()
<< " is not of type "
<< wavePressureFvPatchScalarField::typeName
<< exit(FatalError);
}
const vectorField U(this->U()), Un(this->Un());
const scalarField out(pos0(U & patch().Sf()));
@ -247,32 +242,41 @@ void Foam::waveVelocityFvPatchVectorField::updateCoeffs()
}
else
{
const scalarField& phip =
patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
const vectorField U(this->U());
const scalarField out(pos0(phip));
// Where inflow, fix all velocity components to values specified by the
// wave model.
refValue() = (1 - out)*U;
valueFraction() = (1 - out)*symmTensor::I;
// Where outflow, set the normal component of the velocity to a value
// consistent with phi, but scale it to get the volumentic flow rate
// specified by the wave model. Tangential components are extrapolated.
const scalar QPhip = gSum(out*phip);
const scalar QWave = gSum(out*(U & patch().Sf()));
const vectorField nBySf(patch().Sf()/sqr(patch().magSf()));
if (QPhip > vSmall)
if (inletOutlet_)
{
refValue() += out*(QWave/QPhip)*phip*nBySf;
const scalarField& phip =
patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
const scalarField out(pos0(phip));
// Where inflow, fix all velocity components to values specified by
// the wave model.
refValue() = (1 - out)*U;
valueFraction() = (1 - out)*symmTensor::I;
// Where outflow, set the normal component of the velocity to a
// value consistent with phi, but scale it to get the volumetric
// flow rate specified by the wave model. Tangential components are
// extrapolated.
const scalar QPhip = gSum(out*phip);
const scalar QWave = gSum(out*(U & patch().Sf()));
const vectorField nBySf(patch().Sf()/sqr(patch().magSf()));
if (QPhip > vSmall)
{
refValue() += out*(QWave/QPhip)*phip*nBySf;
}
else
{
refValue() += out*QWave*nBySf;
}
valueFraction() += out*sqr(patch().nf());
}
else
{
refValue() += out*QWave*nBySf;
refValue() = U;
valueFraction() = symmTensor::I;
}
valueFraction() += out*sqr(patch().nf());
}
directionMixedFvPatchVectorField::updateCoeffs();
@ -287,7 +291,8 @@ void Foam::waveVelocityFvPatchVectorField::write
{
directionMixedFvPatchVectorField::write(os);
writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
writeEntryIfDifferent<word>(os, "p", word::null, pName_);
writeEntryIfDifferent<word>(os, "p", "p", pName_);
writeEntryIfDifferent<Switch>(os, "inletOutlet", true, inletOutlet_);
waves_.write(os);
}

47
src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.H
View File

@ -30,32 +30,41 @@ Group
Description
This boundary condition provides a waveVelocity condition. This sets the
velocity to that specified by a superposition of wave models. The
corresponding phase fraction condition looks this condition up and re-uses
the wave modelling.
corresponding phase fraction and pressure conditions look this condition up
and re-use the wave modelling.
Flow reversal will occur in the event that the amplitude of the velocity
oscillation is greater than the mean flow. This triggers special handling,
the form of which depends on whether a pressure field has been specified or
not.
the form of which depends on the inletOutlet flag and whether a wave
pressure condition is being used.
If a pressure field is not specified, then the proportion of the patch over
which the flow is reversed functions in a manner similar to the
flowRateOutletVelocity condition; i.e., the velocity is extrapolated and
then scaled to match the required outlet flow rate. Numerically, this is
still a fixedValue constraint on the normal velocity, just one which tends
to avoid instability. The corresponding pressure condition should be
If a wave pressure condition is not being used, and inletOutlet is false,
then this is a standard fixed value condition, with the value supplied by
the wave model. If flow reversal occurs this state may be unstable. The
corresponding pressure condition should be fixedFluxPressure.
If a wave pressure condition is not being used, and inletOutlet is true or
not specified then the proportion of the patch over which the flow is
reversed functions in a manner similar to the flowRateOutletVelocity
condition; i.e., the velocity is extrapolated and then scaled to match the
required outlet flow rate. Numerically, this is still a fixedValue
constraint on the normal velocity, just one which tends to avoid
instability. Again, the corresponding pressure condition should be
fixedFluxPressure.
If a pressure field is specified, then the normal velocity condition becomes
fixedGradient on outlet faces. This gradient is calculated numerically by
evaluating the wave model on both the patch face and the adjacent cell. The
pressure boundary in this case should be a wavePressure condition. This will
do the opposite; it will fix the pressure value on outlet faces, and the
gradient otherwise.
If a wave pressure condition is being used, then the normal velocity
condition becomes fixedGradient on outlet faces. This gradient is
calculated numerically by evaluating the wave model on both the patch face
and the adjacent cell. The pressure boundary in this case should be a
wavePressure condition. This will do the opposite; it will fix the pressure
value on outlet faces, and the gradient otherwise.
Usage
\table
Property | Description | Req'd? | Default
phi | Name of the flux field | no | phi
p | Name of the pressure field | no | p
inletOutlet | does the condition behave like inletOutlet | no | true
origin | origin of the wave coordinate system | yes |
direction | direction of the mean flow | yes |
speed | speed of the mean flow | yes |
@ -63,8 +72,7 @@ Usage
ramp | ramping function for the mean flow speed | no | None
scale | scale factor along the mean flow direction | no | None
crossScale | scale factor across the mean flow direction | no | None
phi | Name of the flux field | no | phi
p | Name of the pressure field | no |
heightAboveWave | use with the height above the wave | no | false
\endtable
Example of the boundary condition specification:
@ -132,6 +140,9 @@ class waveVelocityFvPatchVectorField
//- Name of the pressure field
const word pName_;
//- Act as an inlet/outlet patch?
const Switch inletOutlet_;
//- Wave superposition
const waveSuperposition waves_;

29
src/waves/waveSuperposition/waveSuperposition.C
View File

@ -182,7 +182,8 @@ Foam::waveSuperposition::waveSuperposition(const objectRegistry& db)
waveAngles_(),
ramp_(),
scale_(),
crossScale_()
crossScale_(),
heightAboveWave_(false)
{}
@ -196,7 +197,8 @@ Foam::waveSuperposition::waveSuperposition(const waveSuperposition& waves)
waveAngles_(waves.waveAngles_),
ramp_(waves.ramp_, false),
scale_(waves.scale_, false),
crossScale_(waves.crossScale_, false)
crossScale_(waves.crossScale_, false),
heightAboveWave_(false)
{}
@ -229,7 +231,8 @@ Foam::waveSuperposition::waveSuperposition
dict.found("crossScale")
? Function1<scalar>::New("crossScale", dict)
: autoPtr<Function1<scalar>>()
)
),
heightAboveWave_(dict.lookupOrDefault<Switch>("heightAboveWave", false))
{
const PtrList<entry> waveEntries(dict.lookup("waves"));
@ -285,6 +288,11 @@ Foam::tmp<Foam::vectorField> Foam::waveSuperposition::ULiquid
vectorField xyz(p.size());
transformation(p, axes, u, xyz);
if (heightAboveWave_)
{
xyz.replace(2, height(t, p));
}
return UMean(t) + (velocity(t, xyz) & axes);
}
@ -302,6 +310,11 @@ Foam::tmp<Foam::vectorField> Foam::waveSuperposition::UGas
axes = tensor(- axes.x(), - axes.y(), axes.z());
if (heightAboveWave_)
{
xyz.replace(2, height(t, p));
}
return UMean(t) + (velocity(t, xyz) & axes);
}
@ -317,6 +330,11 @@ Foam::tmp<Foam::scalarField> Foam::waveSuperposition::pLiquid
vectorField xyz(p.size());
transformation(p, axes, u, xyz);
if (heightAboveWave_)
{
xyz.replace(2, height(t, p));
}
return pressure(t, xyz);
}
@ -358,6 +376,11 @@ void Foam::waveSuperposition::write(Ostream& os) const
{
crossScale_->writeData(os);
}
if (heightAboveWave_)
{
os.writeKeyword("heightAboveWave") << heightAboveWave_
<< token::END_STATEMENT << nl;
}
}

4
src/waves/waveSuperposition/waveSuperposition.H
View File

@ -78,6 +78,10 @@ class waveSuperposition
//- Scaling perpendicular to the flow direction
const autoPtr<Function1<scalar>> crossScale_;
//- Calculate wave properties using the height above the wave (true) or
// the height above the origin (false)?
const Switch heightAboveWave_;
// Private Member Functions