46 const dictionary&
dict,
48 const scalar amplitude,
50 scalar (*modelCelerity)(scalar, scalar, scalar, scalar)
53 const bool haveLength =
dict.found(
"length");
54 const bool havePeriod =
dict.found(
"period");
56 if (haveLength == havePeriod)
59 <<
"Exactly one of either length or period must be specified"
65 return dict.lookup<scalar>(
"length");
69 const scalar period =
dict.lookup<scalar>(
"period");
77 const scalar
length = (length0 + length1)/2;
81 (t < period ? length0 : length1) =
length;
84 return (length0 + length1)/2;
99 return depth*
k(length) >
log(great);
106 const scalar amplitude,
111 return sqrt(g/
k(length)*
tanh(
k(length)*depth));
123 return phase_ +
k()*(
x - celerity()*t);
139 const scalar kzGreat =
log(i*great);
148 const scalar kd =
k()*depth();
162 amplitude_(
wave.amplitude_, false),
163 length_(
wave.length_),
172 const word& modelName,
173 scalar (*modelCelerity)(scalar, scalar, scalar, scalar)
177 depth_(
dict.lookupOrDefault<scalar>(
"depth", great)),
179 length_(length(
dict, depth_, amplitude(), g, modelCelerity)),
180 phase_(
dict.lookup<scalar>(
"phase"))
184 Info<< waveModel::typeName <<
": " << modelName
185 <<
": period = " << length_/
c
186 <<
", length = " << length_ <<
endl;;
204 return amplitude(t)*
cos(angle(t,
x));
214 const scalar ka =
k()*amplitude(t);
Macros for easy insertion into run-time selection tables.
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Run-time selectable general function of one variable.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
A list of keyword definitions, which are a keyword followed by any number of values (e....
A class for managing temporary objects.
Generic base class for waves. Derived classes must implement field functions which return the elevati...
scalar g() const
Get the value of gravity.
virtual void write(Ostream &os) const
Write.
bool deep() const
Return whether shallow and intermediate effects are to be omitted.
virtual ~Airy()
Destructor.
scalar length() const
Get the length.
virtual tmp< scalarField > elevation(const scalar t, const scalarField &x) const
Get the wave elevation at a given time and local coordinates. Local.
scalar k() const
The angular wavenumber [rad/m].
virtual void write(Ostream &os) const
Write.
virtual tmp< vector2DField > velocity(const scalar t, const vector2DField &xz) const
Get the wave velocity at a given time and local coordinates. Local.
scalar amplitude() const
Get the amplitude at steady state.
virtual scalar celerity() const
The wave celerity [m/s].
scalar depth() const
Get the depth.
tmp< scalarField > angle(const scalar t, const scalarField &x) const
Angle of the oscillation [rad].
tmp< vector2DField > vi(const label i, const scalar t, const vector2DField &xz) const
Return the non-dimensionalised i-th harmonic of the velocity.
Airy(const Airy &wave)
Construct a copy.
A class for handling words, derived from string.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
const dimensionedScalar c
Speed of light in a vacuum.
label wave(const fvMesh &mesh, const List< labelPair > &changedPatchAndFaces, const label nCorrections, GeometricField< scalar, PatchField, GeoMesh > &distance, TrackingData &td, GeometricField< DataType, PatchField, GeoMesh > &... data)
Wave distance (and maybe additional) data from faces. If nCorrections is.
defineTypeNameAndDebug(Airy, 0)
addToRunTimeSelectionTable(waveModel, Airy, dictionary)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar exp(const dimensionedScalar &ds)
label log2(label i)
Return the log base 2 by successive bit-shifting of the given label.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar cosh(const dimensionedScalar &ds)
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar tanh(const dimensionedScalar &ds)
dimensionedScalar sinh(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
dimensionedScalar sqrt(const dimensionedScalar &ds)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
tmp< vector2DField > zip(const tmp< scalarField > &x, const tmp< scalarField > &y)
dimensionedScalar cos(const dimensionedScalar &ds)