zhyang-dev/OpenFOAM-5.x
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

setWaves: Corrected handling of multiple wave-type patches

Browse Source
This commit is contained in:
Will Bainbridge
2017-06-02 16:07:36 +01:00
parent 0b7ca1614f
commit b255302fba
7 changed files with 275 additions and 176 deletions
Download Patch File Download Diff File Expand all files Collapse all files

339
applications/utilities/preProcessing/setWaves/setWaves.C
View File

@ -30,111 +30,15 @@ Description
#include "fvCFD.H"
#include "levelSet.H"
#include "pointFields.H"
#include "timeSelector.H"
#include "wallDist.H"
#include "waveAlphaFvPatchScalarField.H"
#include "waveVelocityFvPatchVectorField.H"
#include "waveSuperposition.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
void addWaves
(
const waveSuperposition& waves,
const bool liquid,
volScalarField& alpha,
volVectorField& U
)
{
const scalar t = alpha.time().value();;
const fvMesh& mesh = alpha.mesh();
const pointMesh& pMesh = pointMesh::New(mesh);
// Height fields
const scalarField heightC(waves.height(t, mesh.cellCentres()));
const scalarField heightP(waves.height(t, mesh.points()));
// Velocity fields
const DimensionedField<vector, volMesh>
UGasC
(
IOobject("UGasC", mesh.time().timeName(), mesh),
mesh,
dimVelocity,
waves.UGas(t, mesh.cellCentres())
);
const DimensionedField<vector, pointMesh>
UGasP
(
IOobject("UGasP", mesh.time().timeName(), mesh),
pMesh,
dimVelocity,
waves.UGas(t, mesh.points())
);
const DimensionedField<vector, volMesh>
ULiquidC
(
IOobject("ULiquidC", mesh.time().timeName(), mesh),
mesh,
dimVelocity,
waves.ULiquid(t, mesh.cellCentres())
);
const DimensionedField<vector, pointMesh>
ULiquidP
(
IOobject("ULiquidP", mesh.time().timeName(), mesh),
pMesh,
dimVelocity,
waves.ULiquid(t, mesh.points())
);
// Convert from the level set to volume-averaged fields and sum up
alpha.ref() += levelSetFraction(mesh, heightC, heightP, !liquid);
U.ref() +=
levelSetAverage
(
mesh,
heightC,
heightP,
UGasC,
UGasP,
ULiquidC,
ULiquidP
);
// Now set the boundary fields
forAll(alpha.boundaryField(), patchi)
{
fvPatchScalarField& alphap = alpha.boundaryFieldRef()[patchi];
fvPatchVectorField& Up = U.boundaryFieldRef()[patchi];
const fvPatch& patch = alphap.patch();
// Height fields
const scalarField heightF(waves.height(t, patch.Cf()));
const scalarField heightP(waves.height(t, patch.patch().localPoints()));
// Velocity fields
const vectorField UGasC(waves.UGas(t, mesh.cellCentres()));
const vectorField UGasP(waves.UGas(t, mesh.points()));
const vectorField ULiquidC(waves.ULiquid(t, mesh.cellCentres()));
const vectorField ULiquidP(waves.ULiquid(t, mesh.points()));
alphap == alphap + levelSetFraction(patch, heightF, heightP, !liquid);
Up == Up
+ levelSetAverage
(
patch,
heightC,
heightP,
UGasC,
UGasP,
ULiquidC,
ULiquidP
);
}
}
int main(int argc, char *argv[])
{
timeSelector::addOptions(false, false);
@ -161,6 +65,8 @@ int main(int argc, char *argv[])
#include "createMesh.H"
#include "readGravitationalAcceleration.H"
const pointMesh& pMesh = pointMesh::New(mesh);
forAll(timeDirs, timeI)
{
runTime.setTime(timeDirs[timeI], timeI);
@ -169,7 +75,7 @@ int main(int argc, char *argv[])
mesh.readUpdate();
// Read the phase fraction and velocity fields
// Read the fields which are to be set
volScalarField alpha
(
IOobject
@ -193,20 +99,58 @@ int main(int argc, char *argv[])
mesh
);
// Zero the fields
alpha = dimensionedScalar("0", alpha.dimensions(), 0);
U = dimensionedVector("0", U.dimensions(), vector::zero);
forAll(alpha.boundaryField(), patchi)
{
alpha.boundaryFieldRef()[patchi] == 0;
U.boundaryFieldRef()[patchi] == vector::zero;
}
// Create modelled fields on both cells and points
volScalarField height
(
IOobject("height", runTime.timeName(), mesh),
mesh,
dimensionedScalar("0", dimLength, 0)
);
pointScalarField heightp
(
IOobject("heightp", runTime.timeName(), mesh),
pMesh,
dimensionedScalar("0", dimLength, 0)
);
volVectorField uGas
(
IOobject("uGas", runTime.timeName(), mesh),
mesh,
dimensionedVector("0", dimVelocity, vector::zero)
);
pointVectorField uGasp
(
IOobject("uGasp", runTime.timeName(), mesh),
pMesh,
dimensionedVector("0", dimLength, vector::zero)
);
volVectorField uLiquid
(
IOobject("uLiquid", runTime.timeName(), mesh),
mesh,
dimensionedVector("0", dimVelocity, vector::zero)
);
pointVectorField uLiquidp
(
IOobject("uLiquidp", runTime.timeName(), mesh),
pMesh,
dimensionedVector("0", dimLength, vector::zero)
);
// Loop the patches, looking for wave conditions
forAll(alpha.boundaryField(), patchi)
// The number of wave patches
label nWaves = 0;
// Whether the alpha conditions refer to the liquid phase
bool liquid = false;
// The mean velocity of one of the wave patches
vector UMeanp = vector::zero;
// Loop the patches, averaging and superimposing wave model data
forAll(mesh.boundary(), patchi)
{
const fvPatchScalarField& alphap = alpha.boundaryField()[patchi];
const fvPatchVectorField& Up = U.boundaryField()[patchi];
fvPatchScalarField& alphap = alpha.boundaryFieldRef()[patchi];
fvPatchVectorField& Up = U.boundaryFieldRef()[patchi];
const bool isWave = isA<waveAlphaFvPatchScalarField>(alphap);
@ -223,16 +167,175 @@ int main(int argc, char *argv[])
<< " and vice-versa." << exit(FatalError);
}
if (isWave)
if (!isWave)
{
Info<< "Adding waves from patch " << Up.patch().name() << endl;
addWaves
continue;
}
Info<< "Adding waves from patch " << Up.patch().name() << endl;
const waveSuperposition& waves =
refCast<waveVelocityFvPatchVectorField>(Up).waves();
UMeanp = waves.UMean();
const bool liquidp =
refCast<waveAlphaFvPatchScalarField>(alphap).liquid();
if (nWaves > 0 && liquidp != liquid)
{
FatalErrorInFunction
<< "All " << waveAlphaFvPatchScalarField::typeName
<< "patch fields must be configured for the same phase,"
<< " i.e., the liquid switch must have the same value."
<< exit(FatalError);
}
liquid = liquidp;
const scalar t = runTime.value();
const pointField& ccs = mesh.cellCentres();
const pointField& pts = mesh.points();
// Internal field superposition
height.primitiveFieldRef() += waves.height(t, ccs);
heightp.primitiveFieldRef() += waves.height(t, pts);
uGas.primitiveFieldRef() += waves.UGas(t, ccs) - UMeanp;
uGasp.primitiveFieldRef() += waves.UGas(t, pts) - UMeanp;
uLiquid.primitiveFieldRef() += waves.ULiquid(t, ccs) - UMeanp;
uLiquidp.primitiveFieldRef() += waves.ULiquid(t, pts) - UMeanp;
// Boundary field superposition
forAll(mesh.boundary(), patchj)
{
const pointField& fcs = mesh.boundary()[patchj].Cf();
height.boundaryFieldRef()[patchj] += waves.height(t, fcs);
uGas.boundaryFieldRef()[patchj] += waves.UGas(t, fcs) - UMeanp;
uLiquid.boundaryFieldRef()[patchj] +=
waves.ULiquid(t, fcs) - UMeanp;
}
++ nWaves;
}
// Warn and skip to the next time if no wave patches were found
if (nWaves == 0)
{
WarningInFunction
<< "No " << waveAlphaFvPatchScalarField::typeName << " or "
<< waveVelocityFvPatchVectorField::typeName << " patch fields "
<< "were found. No waves have been set." << endl;
continue;
}
// Create the mean velocity field
volVectorField UMean
(
IOobject("UMean", runTime.timeName(), mesh),
mesh,
dimensionedVector("UMean", dimVelocity, UMeanp)
);
if (nWaves > 1)
{
// Create weighted average fields for the mean velocity
volScalarField weight
(
IOobject("weight", runTime.timeName(), mesh),
mesh,
dimensionedScalar("0", dimless/dimLength, 0)
);
volVectorField weightUMean
(
IOobject("weightUMean", runTime.timeName(), mesh),
mesh,
dimensionedVector("0", dimVelocity/dimLength, vector::zero)
);
// Loop the patches, inverse-distance weighting the mean velocities
forAll(mesh.boundary(), patchi)
{
fvPatchScalarField& alphap = alpha.boundaryFieldRef()[patchi];
fvPatchVectorField& Up = U.boundaryFieldRef()[patchi];
const bool isWave = isA<waveAlphaFvPatchScalarField>(alphap);
if (!isWave)
{
continue;
}
const waveSuperposition& waves =
refCast<waveVelocityFvPatchVectorField>(Up).waves();
UMeanp = waves.UMean();
const volScalarField w
(
refCast<const waveVelocityFvPatchVectorField>(Up).waves(),
refCast<const waveAlphaFvPatchScalarField>(alphap).liquid(),
alpha,
U
1
/(
wallDist(mesh, labelList(1, patchi)).y()
+ dimensionedScalar("ySmall", dimLength, SMALL)
)
);
weight += w;
weightUMean +=
w*dimensionedVector("UMeanp", dimVelocity, UMeanp);
}
// Complete the average for the mean velocity
UMean = weightUMean/weight;
}
// Set the internal fields
alpha.ref() = levelSetFraction(mesh, height, heightp, !liquid);
U.ref() =
UMean
+ levelSetAverage
(
mesh,
height,
heightp,
uGas,
uGasp,
uLiquid,
uLiquidp
);
// Set the boundary fields
forAll(mesh.boundary(), patchi)
{
fvPatchScalarField& alphap = alpha.boundaryFieldRef()[patchi];
fvPatchVectorField& Up = U.boundaryFieldRef()[patchi];
const bool isWave = isA<waveAlphaFvPatchScalarField>(alphap);
if (!isWave)
{
alphap ==
levelSetFraction
(
mesh.boundary()[patchi],
height.boundaryField()[patchi],
heightp.boundaryField()[patchi].patchInternalField(),
!liquid
);
Up ==
UMean.boundaryField()[patchi]
+ levelSetAverage
(
mesh.boundary()[patchi],
height.boundaryField()[patchi],
heightp.boundaryField()[patchi].patchInternalField()(),
uGas.boundaryField()[patchi],
uGasp.boundaryField()[patchi].patchInternalField()(),
uLiquid.boundaryField()[patchi],
uLiquidp.boundaryField()[patchi].patchInternalField()()
);
}
else
{
alphap == refCast<waveAlphaFvPatchScalarField>(alphap).alpha();
Up == refCast<waveVelocityFvPatchVectorField>(Up).U();
}
}

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

@ -30,30 +30,6 @@ License
#include "surfaceFields.H"
#include "volFields.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::tmp<Foam::scalarField> Foam::waveAlphaFvPatchScalarField::alpha() const
{
const scalar t = db().time().timeOutputValue();
const waveVelocityFvPatchVectorField& Up =
refCast<const waveVelocityFvPatchVectorField>
(
patch().lookupPatchField<volVectorField, scalar>(UName_)
);
const waveSuperposition& waves = Up.waves();
return
levelSetFraction
(
patch(),
waves.height(t, patch().Cf()),
waves.height(t, patch().patch().localPoints()),
!liquid_
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::waveAlphaFvPatchScalarField::waveAlphaFvPatchScalarField
@ -141,6 +117,28 @@ Foam::waveAlphaFvPatchScalarField::waveAlphaFvPatchScalarField
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::tmp<Foam::scalarField> Foam::waveAlphaFvPatchScalarField::alpha() const
{
const scalar t = db().time().timeOutputValue();
const waveVelocityFvPatchVectorField& Up =
refCast<const waveVelocityFvPatchVectorField>
(
patch().lookupPatchField<volVectorField, scalar>(UName_)
);
const waveSuperposition& waves = Up.waves();
return
levelSetFraction
(
patch(),
waves.height(t, patch().Cf()),
waves.height(t, patch().patch().localPoints()),
!liquid_
);
}
void Foam::waveAlphaFvPatchScalarField::updateCoeffs()
{
if (updated())

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

@ -89,12 +89,6 @@ class waveAlphaFvPatchScalarField
const bool inletOutlet_;
// Private Member Functions
//- Return the current modelled phase fraction field
tmp<scalarField> alpha() const;
public:
//- Runtime type information
@ -175,6 +169,9 @@ public:
// Evaluation functions
//- Return the current modelled phase fraction field
tmp<scalarField> alpha() const;
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();

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

@ -29,26 +29,6 @@ License
#include "surfaceFields.H"
#include "volFields.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::tmp<Foam::vectorField> Foam::waveVelocityFvPatchVectorField::U() const
{
const scalar t = db().time().timeOutputValue();
return
levelSetAverage
(
patch(),
waves_.height(t, patch().Cf()),
waves_.height(t, patch().patch().localPoints()),
waves_.UGas(t, patch().Cf())(),
waves_.UGas(t, patch().patch().localPoints())(),
waves_.ULiquid(t, patch().Cf())(),
waves_.ULiquid(t, patch().patch().localPoints())()
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField
@ -127,6 +107,24 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::tmp<Foam::vectorField> Foam::waveVelocityFvPatchVectorField::U() const
{
const scalar t = db().time().timeOutputValue();
return
levelSetAverage
(
patch(),
waves_.height(t, patch().Cf()),
waves_.height(t, patch().patch().localPoints()),
waves_.UGas(t, patch().Cf())(),
waves_.UGas(t, patch().patch().localPoints())(),
waves_.ULiquid(t, patch().Cf())(),
waves_.ULiquid(t, patch().patch().localPoints())()
);
}
void Foam::waveVelocityFvPatchVectorField::updateCoeffs()
{
if (updated())

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

@ -104,12 +104,6 @@ class waveVelocityFvPatchVectorField
const waveSuperposition waves_;
// Private Member Functions
//- Return the current modelled velocity field
tmp<vectorField> U() const;
public:
//- Runtime type information
@ -190,6 +184,9 @@ public:
// Evaluation functions
//- Return the current modelled velocity field
tmp<vectorField> U() const;
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();

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

@ -245,7 +245,7 @@ Foam::tmp<Foam::vectorField> Foam::waveSuperposition::ULiquid
vectorField xyz(p.size());
transformation(p, axes, u, xyz);
return direction_*speed_ + (velocity(t, xyz) & axes);
return UMean() + (velocity(t, xyz) & axes);
}
@ -262,7 +262,7 @@ Foam::tmp<Foam::vectorField> Foam::waveSuperposition::UGas
axes = tensor(- axes.x(), - axes.y(), axes.z());
return direction_*speed_ + (velocity(t, xyz) & axes);
return UMean() + (velocity(t, xyz) & axes);
}

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

@ -139,6 +139,12 @@ public:
//- Get the gas velocity at a given time and global positions
tmp<vectorField> UGas(const scalar t, const vectorField& p) const;
//- Get the mean flow velocity
inline vector UMean() const
{
return direction_*speed_;
}
//- Write
void write(Ostream&) const;
};