From 5925868fb7a54f3023d2ccc962079bfad03c9972 Mon Sep 17 00:00:00 2001 From: Will Bainbridge Date: Tue, 18 Dec 2018 10:20:28 +0000 Subject: [PATCH] waves: Moved mean velocity specification back into the wave models With the inclusion of boundary layer modelling in the gas, the separation of wave perturbation from and mean flow became less useful, and potentially prevents further extension to support similar boundary layer modelling in the liquid. The mean velocity entry, UMean, is now needed in the constant/waveProperties file rather than in the waveVelocity boundary condition. --- .../preProcessing/setWaves/setWaves.C | 29 ++++++--------- etc/caseDicts/annotated/setWavesDict | 3 -- .../waveAlpha/waveAlphaFvPatchScalarField.C | 13 +++---- .../wavePressureFvPatchScalarField.C | 29 +++++++-------- .../waveVelocityFvPatchVectorField.C | 36 ++++++++----------- .../waveVelocityFvPatchVectorField.H | 17 --------- .../waveAtmBoundaryLayerSuperposition.C | 33 +++++++++++------ .../waveAtmBoundaryLayerSuperposition.H | 27 +++++++------- .../waveSuperposition/waveSuperposition.C | 19 +++++----- .../waveSuperposition/waveSuperposition.H | 26 ++++++-------- .../interFoam/RAS/DTCHullWave/0/U.orig | 7 ++-- .../RAS/DTCHullWave/constant/waveProperties | 4 +++ .../RAS/DTCHullWave/system/setWavesDict | 7 ---- .../interFoam/laminar/wave/0/U.orig | 1 - .../laminar/wave/constant/waveProperties | 2 ++ .../laminar/wave/system/blockMeshDict | 2 +- .../laminar/wave/system/setWavesDict | 2 -- 17 files changed, 106 insertions(+), 151 deletions(-) diff --git a/applications/utilities/preProcessing/setWaves/setWaves.C b/applications/utilities/preProcessing/setWaves/setWaves.C index 85438a2496..bc48d78bbf 100644 --- a/applications/utilities/preProcessing/setWaves/setWaves.C +++ b/applications/utilities/preProcessing/setWaves/setWaves.C @@ -91,12 +91,6 @@ int main(int argc, char *argv[]) const word alphaName = setWavesDict.lookupOrDefault("alpha", "alpha"); const word UName = setWavesDict.lookupOrDefault("U", "U"); const bool liquid = setWavesDict.lookupOrDefault("liquid", true); - const dimensionedVector UMean - ( - "UMean", - dimVelocity, - setWavesDict.lookup("UMean") - ); // Get the wave models const waveSuperposition& waves = waveSuperposition::New(mesh); @@ -172,33 +166,30 @@ int main(int argc, char *argv[]) dimensionedVector("0", dimVelocity, vector::zero) ); - // Offset - const vector offset = UMean.value()*t; - // Cell centres and points const pointField& ccs = mesh.cellCentres(); const pointField& pts = mesh.points(); // Internal field - h.primitiveFieldRef() = waves.height(t, ccs + offset); - hp.primitiveFieldRef() = waves.height(t, pts + offset); - uGas.primitiveFieldRef() = waves.UGas(t, ccs + offset); - uGasp.primitiveFieldRef() = waves.UGas(t, pts + offset); - uLiq.primitiveFieldRef() = waves.ULiquid(t, ccs + offset); - uLiqp.primitiveFieldRef() = waves.ULiquid(t, pts + offset); + h.primitiveFieldRef() = waves.height(t, ccs); + hp.primitiveFieldRef() = waves.height(t, pts); + uGas.primitiveFieldRef() = waves.UGas(t, ccs); + uGasp.primitiveFieldRef() = waves.UGas(t, pts); + uLiq.primitiveFieldRef() = waves.ULiquid(t, ccs); + uLiqp.primitiveFieldRef() = waves.ULiquid(t, pts); // Boundary fields forAll(mesh.boundary(), patchj) { const pointField& fcs = mesh.boundary()[patchj].Cf(); - h.boundaryFieldRef()[patchj] = waves.height(t, fcs + offset); - uGas.boundaryFieldRef()[patchj] = waves.UGas(t, fcs + offset); - uLiq.boundaryFieldRef()[patchj] = waves.ULiquid(t, fcs + offset); + h.boundaryFieldRef()[patchj] = waves.height(t, fcs); + uGas.boundaryFieldRef()[patchj] = waves.UGas(t, fcs); + uLiq.boundaryFieldRef()[patchj] = waves.ULiquid(t, fcs); } // Set the fields alpha == levelSetFraction(h, hp, !liquid); - U == UMean + levelSetAverage(h, hp, uGas, uGasp, uLiq, uLiqp); + U == levelSetAverage(h, hp, uGas, uGasp, uLiq, uLiqp); // Set the boundary fields forAll(mesh.boundary(), patchi) diff --git a/etc/caseDicts/annotated/setWavesDict b/etc/caseDicts/annotated/setWavesDict index b603c25136..9763349781 100644 --- a/etc/caseDicts/annotated/setWavesDict +++ b/etc/caseDicts/annotated/setWavesDict @@ -24,8 +24,5 @@ U U; // under the waves (true) or the gas over the waves (false) liquid true; -// The mean flow velocity over which to superimpose waves -UMean (1 0 0); - // ************************************************************************* // diff --git a/src/waves/derivedFvPatchFields/waveAlpha/waveAlphaFvPatchScalarField.C b/src/waves/derivedFvPatchFields/waveAlpha/waveAlphaFvPatchScalarField.C index c1ee40ee97..d0ca44dbb6 100644 --- a/src/waves/derivedFvPatchFields/waveAlpha/waveAlphaFvPatchScalarField.C +++ b/src/waves/derivedFvPatchFields/waveAlpha/waveAlphaFvPatchScalarField.C @@ -123,18 +123,13 @@ Foam::tmp Foam::waveAlphaFvPatchScalarField::alpha() const { const scalar t = db().time().timeOutputValue(); const waveSuperposition& waves = waveSuperposition::New(db()); - const waveVelocityFvPatchVectorField& Up = - refCast - ( - patch().lookupPatchField(UName_) - ); return levelSetFraction ( patch(), - waves.height(t, patch().Cf() + Up.offset()), - waves.height(t, patch().patch().localPoints() + Up.offset()), + waves.height(t, patch().Cf()), + waves.height(t, patch().patch().localPoints()), !liquid_ ); } @@ -159,8 +154,8 @@ Foam::tmp Foam::waveAlphaFvPatchScalarField::alphan() const levelSetFraction ( meshs, - waves.height(t, meshs.cellCentres() + Up.offset())(), - waves.height(t, meshs.points() + Up.offset())(), + waves.height(t, meshs.cellCentres())(), + waves.height(t, meshs.points())(), !liquid_ ) ); diff --git a/src/waves/derivedFvPatchFields/wavePressure/wavePressureFvPatchScalarField.C b/src/waves/derivedFvPatchFields/wavePressure/wavePressureFvPatchScalarField.C index 7fed7ae8f9..c89074266d 100644 --- a/src/waves/derivedFvPatchFields/wavePressure/wavePressureFvPatchScalarField.C +++ b/src/waves/derivedFvPatchFields/wavePressure/wavePressureFvPatchScalarField.C @@ -117,22 +117,17 @@ Foam::tmp Foam::wavePressureFvPatchScalarField::p() const { const scalar t = db().time().timeOutputValue(); const waveSuperposition& waves = waveSuperposition::New(db()); - const waveVelocityFvPatchVectorField& Up = - refCast - ( - patch().lookupPatchField(UName_) - ); return levelSetAverage ( patch(), - waves.height(t, patch().Cf() + Up.offset()), - waves.height(t, patch().patch().localPoints() + Up.offset()), - waves.pGas(t, patch().Cf() + Up.offset())(), - waves.pGas(t, patch().patch().localPoints() + Up.offset())(), - waves.pLiquid(t, patch().Cf() + Up.offset())(), - waves.pLiquid(t, patch().patch().localPoints() + Up.offset())() + waves.height(t, patch().Cf()), + waves.height(t, patch().patch().localPoints()), + waves.pGas(t, patch().Cf())(), + waves.pGas(t, patch().patch().localPoints())(), + waves.pLiquid(t, patch().Cf())(), + waves.pLiquid(t, patch().patch().localPoints())() ); } @@ -156,12 +151,12 @@ Foam::tmp Foam::wavePressureFvPatchScalarField::pn() const levelSetAverage ( meshs, - waves.height(t, meshs.cellCentres() + Up.offset())(), - waves.height(t, meshs.points() + Up.offset())(), - waves.pGas(t, meshs.cellCentres() + Up.offset())(), - waves.pGas(t, meshs.points() + Up.offset())(), - waves.pLiquid(t, meshs.cellCentres() + Up.offset())(), - waves.pLiquid(t, meshs.points() + Up.offset())() + waves.height(t, meshs.cellCentres())(), + waves.height(t, meshs.points())(), + waves.pGas(t, meshs.cellCentres())(), + waves.pGas(t, meshs.points())(), + waves.pLiquid(t, meshs.cellCentres())(), + waves.pLiquid(t, meshs.points())() ) ); diff --git a/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.C b/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.C index 75a34f4b84..c22c180d79 100644 --- a/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.C +++ b/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.C @@ -42,7 +42,6 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField phiName_("phi"), pName_("p"), inletOutlet_(true), - UMean_(nullptr), faceCellSubset_(nullptr), faceCellSubsetTimeIndex_(-1) { @@ -63,7 +62,6 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField phiName_(dict.lookupOrDefault("phi", "phi")), pName_(dict.lookupOrDefault("p", "p")), inletOutlet_(dict.lookupOrDefault("inletOutlet", true)), - UMean_(Function1::New("UMean", dict)), faceCellSubset_(nullptr), faceCellSubsetTimeIndex_(-1) { @@ -94,7 +92,6 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField phiName_(ptf.phiName_), pName_(ptf.pName_), inletOutlet_(ptf.inletOutlet_), - UMean_(ptf.UMean_, false), faceCellSubset_(nullptr), faceCellSubsetTimeIndex_(-1) {} @@ -109,7 +106,6 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField phiName_(ptf.phiName_), pName_(ptf.pName_), inletOutlet_(ptf.inletOutlet_), - UMean_(ptf.UMean_, false), faceCellSubset_(nullptr), faceCellSubsetTimeIndex_(-1) {} @@ -125,7 +121,6 @@ Foam::waveVelocityFvPatchVectorField::waveVelocityFvPatchVectorField phiName_(ptf.phiName_), pName_(ptf.pName_), inletOutlet_(ptf.inletOutlet_), - UMean_(ptf.UMean_, false), faceCellSubset_(nullptr), faceCellSubsetTimeIndex_(-1) {} @@ -165,16 +160,15 @@ Foam::tmp Foam::waveVelocityFvPatchVectorField::U() const const waveSuperposition& waves = waveSuperposition::New(db()); return - UMean_->value(t) - + levelSetAverage + levelSetAverage ( patch(), - waves.height(t, patch().Cf() + offset()), - waves.height(t, patch().patch().localPoints() + offset()), - waves.UGas(t, patch().Cf() + offset())(), - waves.UGas(t, patch().patch().localPoints() + offset())(), - waves.ULiquid(t, patch().Cf() + offset())(), - waves.ULiquid(t, patch().patch().localPoints() + offset())() + 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())() ); } @@ -190,16 +184,15 @@ Foam::tmp Foam::waveVelocityFvPatchVectorField::Un() const const vectorField Us ( - UMean_->value(t) - + levelSetAverage + levelSetAverage ( meshs, - waves.height(t, meshs.cellCentres() + offset())(), - waves.height(t, meshs.points() + offset())(), - waves.UGas(t, meshs.cellCentres() + offset())(), - waves.UGas(t, meshs.points() + offset())(), - waves.ULiquid(t, meshs.cellCentres() + offset())(), - waves.ULiquid(t, meshs.points() + offset())() + waves.height(t, meshs.cellCentres())(), + waves.height(t, meshs.points())(), + waves.UGas(t, meshs.cellCentres())(), + waves.UGas(t, meshs.points())(), + waves.ULiquid(t, meshs.cellCentres())(), + waves.ULiquid(t, meshs.points())() ) ); @@ -297,7 +290,6 @@ void Foam::waveVelocityFvPatchVectorField::write writeEntryIfDifferent(os, "phi", "phi", phiName_); writeEntryIfDifferent(os, "p", "p", pName_); writeEntryIfDifferent(os, "inletOutlet", true, inletOutlet_); - UMean_->writeData(os); } diff --git a/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.H b/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.H index 71de15223a..73e2be476f 100644 --- a/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.H +++ b/src/waves/derivedFvPatchFields/waveVelocity/waveVelocityFvPatchVectorField.H @@ -62,7 +62,6 @@ Usage 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 - UMean | velocity of the mean flow | yes | ramp | ramping function for the mean flow speed | no | None \endtable @@ -74,7 +73,6 @@ Usage phi phi; p p; inletOutlet yes; - UMean (2 0 0); ramp constant 1; } \endverbatim @@ -119,9 +117,6 @@ class waveVelocityFvPatchVectorField //- Act as an inlet/outlet patch? const Switch inletOutlet_; - //- Mean velocity - const autoPtr> UMean_; - //- Mesh subset corresponding to the patch adjacent cells mutable autoPtr faceCellSubset_; @@ -218,18 +213,6 @@ public: //- Return the current modelled velocity field in the neighbour cell tmp Un() const; - //- Return the mean velocity - inline vector UMean() const - { - return UMean_->value(db().time().timeOutputValue()); - } - - //- Return the distance offset for the wave models - inline vector offset() const - { - return - UMean_->integrate(0, db().time().timeOutputValue()); - } - //- Update the coefficients associated with the patch field virtual void updateCoeffs(); diff --git a/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.C b/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.C index eb07ed400d..8b79513a51 100644 --- a/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.C +++ b/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.C @@ -77,19 +77,30 @@ Foam::tmp Foam::waveAtmBoundaryLayerSuperposition::UGas db().lookupObject("g").value() ); - const scalar h0 = - gHat & origin(); + const scalar h0 = - gHat & origin_; - atmBoundaryLayer atm - ( - normalised(UGasRef_), - - gHat, - mag(UGasRef_), - h0 + hRef_, - scalarField(p.size(), hWaveMax_ - hWaveMin_), - scalarField(p.size(), h0 + hWaveMin_) - ); + const vector UGasRefRel = UGasRef_ - UMean_->value(t); - return waveSuperposition::UGas(t, p) + atm.U(p); + const scalar magUGasRefRel = mag(UGasRefRel); + + tmp tU = waveSuperposition::UGas(t, p); + + if (magUGasRefRel > 0) + { + atmBoundaryLayer atm + ( + UGasRefRel/magUGasRefRel, + - gHat, + magUGasRefRel, + h0 + hRef_, + scalarField(p.size(), hWaveMax_ - hWaveMin_), + scalarField(p.size(), h0 + hWaveMin_) + ); + + tU.ref() += atm.U(p); + } + + return tU; } diff --git a/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.H b/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.H index 4e3b443a04..47d68047ba 100644 --- a/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.H +++ b/src/waves/waveSuperpositions/waveAtmBoundaryLayerSuperposition/waveAtmBoundaryLayerSuperposition.H @@ -26,24 +26,23 @@ Class Description An extension of waveSuperposition which adds an atmospheric boundary layer - model to the gas velocity. The user supplies a gas velocity relative to the - mean liquid velocity and a height above the wave coordinate system origin - at which that velocity is reached. Also needed are a maximum and minimum - wave height which are used to set the surface roughness in the boundary - layer model. It is not trivial to determine these from an arbitrary - superposition of differently oriented wave models, so they are required as - user inputs instead. For a pure sinusoidal wave, the maximum and minimum - wave heights can be set to positive and negative amplitude, respectively. + model to the gas velocity. The user supplies a gas velocity and a height + above the wave coordinate system origin at which that velocity is reached. + Also needed are a maximum and minimum wave height which are used to set the + surface roughness in the boundary layer model. It is not trivial to + determine these from an arbitrary superposition of differently oriented + wave models, so they are required as user inputs instead. For a pure + sinusoidal wave, the maximum and minimum wave heights can be set to + positive and negative amplitude, respectively. Usage \table - Property | Description | Req'd? | Default - UGasRef | The gas velocity relative to the \\ - liquid at the reference height | yes | + Property | Description | Req'd? | Default + UGasRef | The gas velocity at the reference height | yes | hRef | The reference height relative to the \\ - origin of the wave coordinate system | yes | - hWaveMin | The minimum wave elevation | yes | - hWaveMax | The maximum wave elevation | yes | + origin of the wave coordinate system | yes | + hWaveMin | The minimum wave elevation | yes | + hWaveMax | The maximum wave elevation | yes | \endtable Example specification: diff --git a/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.C b/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.C index 19489e7b6f..d8d9caf921 100644 --- a/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.C +++ b/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.C @@ -48,6 +48,7 @@ namespace Foam void Foam::waveSuperposition::transformation ( + const scalar t, const vectorField& p, tensor& axes, vectorField& xyz @@ -64,7 +65,7 @@ void Foam::waveSuperposition::transformation axes = tensor(dSurfHat, - gHat ^ dSurfHat, - gHat); - xyz = axes & (p - origin_); + xyz = axes & (p - origin_ - UMean_->integrate(0, t)); } @@ -204,6 +205,7 @@ Foam::waveSuperposition::waveSuperposition(const objectRegistry& db) direction_(lookup("direction")), waveModels_(), waveAngles_(), + UMean_(Function1::New("UMean", *this)), scale_ ( found("scale") @@ -252,7 +254,7 @@ Foam::tmp Foam::waveSuperposition::height { tensor axes; vectorField xyz(p.size()); - transformation(p, axes, xyz); + transformation(t, p, axes, xyz); return xyz.component(2) @@ -268,14 +270,14 @@ Foam::tmp Foam::waveSuperposition::ULiquid { tensor axes; vectorField xyz(p.size()); - transformation(p, axes, xyz); + transformation(t, p, axes, xyz); if (heightAboveWave_) { xyz.replace(2, height(t, p)); } - return velocity(t, xyz) & axes; + return UMean_->value(t) + (velocity(t, xyz) & axes); } @@ -287,7 +289,7 @@ Foam::tmp Foam::waveSuperposition::UGas { tensor axes; vectorField xyz(p.size()); - transformation(p, axes, xyz); + transformation(t, p, axes, xyz); axes = tensor(- axes.x(), - axes.y(), axes.z()); @@ -298,7 +300,7 @@ Foam::tmp Foam::waveSuperposition::UGas xyz.replace(2, - xyz.component(2)); - return velocity(t, xyz) & axes; + return UMean_->value(t) + (velocity(t, xyz) & axes); } @@ -310,7 +312,7 @@ Foam::tmp Foam::waveSuperposition::pLiquid { tensor axes; vectorField xyz(p.size()); - transformation(p, axes, xyz); + transformation(t, p, axes, xyz); if (heightAboveWave_) { @@ -329,7 +331,7 @@ Foam::tmp Foam::waveSuperposition::pGas { tensor axes; vectorField xyz(p.size()); - transformation(p, axes, xyz); + transformation(t, p, axes, xyz); axes = tensor(- axes.x(), - axes.y(), axes.z()); @@ -358,6 +360,7 @@ void Foam::waveSuperposition::write(Ostream& os) const << nl << decrIndent << indent << token::END_BLOCK << nl; } os << decrIndent << token::END_LIST << token::END_STATEMENT << nl; + UMean_->writeData(os); if (scale_.valid()) { scale_->writeData(os); diff --git a/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.H b/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.H index a332d96aa8..d001a8475d 100644 --- a/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.H +++ b/src/waves/waveSuperpositions/waveSuperposition/waveSuperposition.H @@ -36,6 +36,7 @@ Usage origin | origin of the wave coordinate system | yes | direction | direction of the wave coordinate system | yes | waves | list of wave models to superimpose | yes | + UMean | velocity of the mean flow | yes | scale | scale factor in the direction | no | None crossScale | scale factor perpendicular to the direction | no | None heightAboveWave | use the height above the wave as the vertical \\ @@ -63,6 +64,7 @@ Usage angle 0; } ); + UMean (2 0 0); scale table ((100 1) (200 0)); crossScale constant 1; heightAboveWave no; @@ -92,7 +94,9 @@ class waveSuperposition : public IOdictionary { - // Private Data +protected: + + // Protected Data //- The origin of the wave coordinate system const vector origin_; @@ -106,6 +110,9 @@ class waveSuperposition //- The angle relative to the direction at which the waves propagate scalarList waveAngles_; + //- Mean velocity + const autoPtr> UMean_; + //- Scaling in the local x-direction const autoPtr> scale_; @@ -117,11 +124,12 @@ class waveSuperposition const Switch heightAboveWave_; - // Private Member Functions + // Protected Member Functions //- Get the transformation to actual coordinates void transformation ( + const scalar t, const vectorField& p, tensor& axes, vectorField& xyz @@ -193,20 +201,6 @@ public: // Member Functions - // Access - - //- Return the origin of the wave coordinate system - const vector& origin() const - { - return origin_; - } - - //- Return the direction of the wave coordinate system - const vector& direction() const - { - return direction_; - } - //- Get the height above the waves at a given time and global positions virtual tmp height ( diff --git a/tutorials/multiphase/interFoam/RAS/DTCHullWave/0/U.orig b/tutorials/multiphase/interFoam/RAS/DTCHullWave/0/U.orig index f0e3d4678b..69aeb8b99d 100644 --- a/tutorials/multiphase/interFoam/RAS/DTCHullWave/0/U.orig +++ b/tutorials/multiphase/interFoam/RAS/DTCHullWave/0/U.orig @@ -15,11 +15,11 @@ FoamFile } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -UMean 1.668; +#include "$FOAM_CASE/constant/waveProperties"; dimensions [0 1 -1 0 0 0 0]; -internalField uniform (#neg $UMean 0 0); +internalField uniform ($UxMean 0 0); boundaryField { @@ -29,14 +29,13 @@ boundaryField inlet { type waveVelocity; - UMean (#neg $UMean 0 0); } outlet { type outletPhaseMeanVelocity; alpha alpha.water; - UnMean $UMean; + UnMean #neg $UxMean; value $internalField; } diff --git a/tutorials/multiphase/interFoam/RAS/DTCHullWave/constant/waveProperties b/tutorials/multiphase/interFoam/RAS/DTCHullWave/constant/waveProperties index 38c13b9aae..8fbf6c89ee 100644 --- a/tutorials/multiphase/interFoam/RAS/DTCHullWave/constant/waveProperties +++ b/tutorials/multiphase/interFoam/RAS/DTCHullWave/constant/waveProperties @@ -30,6 +30,10 @@ waves } ); +UxMean -1.668; + +UMean ($UxMean 0 0); + scale table ((4 1) (12 0)); diff --git a/tutorials/multiphase/interFoam/RAS/DTCHullWave/system/setWavesDict b/tutorials/multiphase/interFoam/RAS/DTCHullWave/system/setWavesDict index 362f80e95f..969ef86bf3 100644 --- a/tutorials/multiphase/interFoam/RAS/DTCHullWave/system/setWavesDict +++ b/tutorials/multiphase/interFoam/RAS/DTCHullWave/system/setWavesDict @@ -14,14 +14,7 @@ FoamFile } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -U.orig -{ - #include "$FOAM_CASE/0/U.orig"; -} - alpha alpha.water; -UMean $U.orig.boundaryField.inlet.UMean; - // ************************************************************************* // diff --git a/tutorials/multiphase/interFoam/laminar/wave/0/U.orig b/tutorials/multiphase/interFoam/laminar/wave/0/U.orig index a5ca3c0971..4a7fe020eb 100644 --- a/tutorials/multiphase/interFoam/laminar/wave/0/U.orig +++ b/tutorials/multiphase/interFoam/laminar/wave/0/U.orig @@ -26,7 +26,6 @@ boundaryField left { type waveVelocity; - UMean (2 0 0); } right { diff --git a/tutorials/multiphase/interFoam/laminar/wave/constant/waveProperties b/tutorials/multiphase/interFoam/laminar/wave/constant/waveProperties index 9daba32311..9f73c7b556 100644 --- a/tutorials/multiphase/interFoam/laminar/wave/constant/waveProperties +++ b/tutorials/multiphase/interFoam/laminar/wave/constant/waveProperties @@ -30,6 +30,8 @@ waves } ); +UMean (2 0 0); + scale table ((1200 1) (1800 0)); crossScale constant 1; diff --git a/tutorials/multiphase/interFoam/laminar/wave/system/blockMeshDict b/tutorials/multiphase/interFoam/laminar/wave/system/blockMeshDict index 5c7ba6693a..6f0b746720 100644 --- a/tutorials/multiphase/interFoam/laminar/wave/system/blockMeshDict +++ b/tutorials/multiphase/interFoam/laminar/wave/system/blockMeshDict @@ -65,7 +65,7 @@ boundary } bottom { - type patch; + type wall; faces ( (0 1 5 4) diff --git a/tutorials/multiphase/interFoam/laminar/wave/system/setWavesDict b/tutorials/multiphase/interFoam/laminar/wave/system/setWavesDict index 8d9d7a085b..969ef86bf3 100644 --- a/tutorials/multiphase/interFoam/laminar/wave/system/setWavesDict +++ b/tutorials/multiphase/interFoam/laminar/wave/system/setWavesDict @@ -16,7 +16,5 @@ FoamFile alpha alpha.water; -UMean (2 0 0); - // ************************************************************************* //