48 int main(
int argc,
char *argv[])
59 "name of the volume fraction field, default is \"alpha\""
66 "name of the velocity field, default is \"U\""
72 "the volume fraction field is that of the gas above the wave"
84 #include "readGravitationalAcceleration.H"
89 const word alphaName = setWavesDict.lookupOrDefault<
word>(
"alpha",
"alpha");
90 const word UName = setWavesDict.lookupOrDefault<
word>(
"U",
"U");
91 const bool liquid = setWavesDict.lookupOrDefault<
bool>(
"liquid",
true);
98 runTime.setTime(
timeDirs[timeI], timeI);
99 const scalar t = runTime.value();
101 Info<<
"Time = " << runTime.userTimeName() <<
nl <<
endl;
132 IOobject(
"h", runTime.name(), mesh),
138 IOobject(
"hp", runTime.name(), mesh),
144 IOobject(
"uGas", runTime.name(), mesh),
150 IOobject(
"uGasp", runTime.name(), mesh),
156 IOobject(
"uLiq", runTime.name(), mesh),
162 IOobject(
"uLiqp", runTime.name(), mesh),
172 h.primitiveFieldRef() = waves.height(t, ccs);
173 hp.primitiveFieldRef() = waves.height(t, pts);
174 uGas.primitiveFieldRef() = waves.UGas(t, ccs);
175 uGasp.primitiveFieldRef() = waves.UGas(t, pts);
176 uLiq.primitiveFieldRef() = waves.ULiquid(t, ccs);
177 uLiqp.primitiveFieldRef() = waves.ULiquid(t, pts);
180 forAll(mesh.boundary(), patchj)
182 const pointField& fcs = mesh.boundary()[patchj].Cf();
183 h.boundaryFieldRef()[patchj] = waves.height(t, fcs);
184 uGas.boundaryFieldRef()[patchj] = waves.UGas(t, fcs);
185 uLiq.boundaryFieldRef()[patchj] = waves.ULiquid(t, fcs);
201 !isA<wallPolyPatch>(patch)
202 || isA<waveAlphaFvPatchScalarField>(alphap)
203 || isA<waveVelocityFvPatchVectorField>(Up)
206 alphap == alphaNoBCs.boundaryField()[
patchi];
207 Up == UNoBCs.boundaryField()[
patchi];
#define forAll(list, i)
Loop across all elements in list.
static pointMesh & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
Generic GeometricField class.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
const word & name() const
Return name.
static void addOption(const word &opt, const string ¶m="", const string &usage="")
Add to an option to validOptions with usage information.
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Mesh representing a set of points created from polyMesh.
A patch is a list of labels that address the faces in the global face list.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
virtual bool write(const bool write=true) const
Write using setting from DB.
static void addOptions(const bool constant=true, const bool withZero=false)
Add the options handled by timeSelector to argList::validOptions.
static instantList selectIfPresent(Time &runTime, const argList &args)
If any time option provided return the set of times (as select0)
A wrapper around a list of wave models. Superimposes the modelled values of elevation and velocity....
static const waveSuperposition & New(const objectRegistry &db)
Return a reference to the wave model on the given database,.
A class for handling words, derived from string.
int main(int argc, char *argv[])
static instantList timeDirs
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const dimensionedScalar h
Planck constant.
tmp< scalarField > levelSetFraction(const fvMesh &mesh, const scalarField &levelC, const scalarField &levelP, const bool above)
Calculate the volume-fraction that a level set occupies. This gives the.
IOdictionary systemDict(const word &dictName, const argList &args, const objectRegistry &ob, const word ®ionName=polyMesh::defaultRegion)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const dimensionSet dimLength
const dimensionSet dimVelocity
tmp< Field< Type > > levelSetAverage(const fvMesh &mesh, const scalarField &levelC, const scalarField &levelP, const Field< Type > &positiveC, const Field< Type > &positiveP, const Field< Type > &negativeC, const Field< Type > &negativeP)
Calculate the average value of two fields, one on each side of a level set.
Foam::argList args(argc, argv)