46 const scalar amplitude,
51 static const scalar kdGreat =
log(great);
58 const scalar C2ByC0 = (2 + 7*
sqr(
S))/4/
sqr(1 -
S);
75 const word& modelName,
76 scalar (*modelCelerity)(scalar, scalar, scalar, scalar)
79 Airy(
dict, g, modelName, modelCelerity)
97 static const scalar kdGreat =
log(great);
98 const scalar kd =
min(
max(
k()*depth(), - kdGreat), kdGreat);
99 const scalar ka =
k()*amplitude(t);
101 const scalar
T = deep() ? 1 :
tanh(kd);
103 const scalar B22 = (3/
sqr(
T) - 1)/
T/4;
112 + (1/
k())*
sqr(ka)*B22*
cos(2*angle(t,
x));
122 static const scalar kdGreat =
log(great);
123 const scalar kd =
min(
max(
k()*depth(), - kdGreat), kdGreat);
124 const scalar ka =
k()*amplitude(t);
126 const scalar A22ByA11 = deep() ? 0 : 0.375/
pow3(
sinh(kd));
130 const scalar A11 = 1/
sinh(kd);
131 Info<<
"A22 = " << A22ByA11*A11 <<
endl;
Macros for easy insertion into run-time selection tables.
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.
bool deep() const
Return whether shallow and intermediate effects are to be omitted.
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 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.
virtual ~Stokes2()
Destructor.
Stokes2(const dictionary &dict, const scalar g, const word &modelName=Stokes2::typeName, scalar(*modelCelerity)(scalar, scalar, scalar, scalar)=&Stokes2::celerity)
Construct from a dictionary and gravity.
virtual tmp< scalarField > elevation(const scalar t, const scalarField &x) const
Get the wave elevation at a given time and local coordinates. Local.
virtual tmp< vector2DField > velocity(const scalar t, const vector2DField &xz) const
Get the wave velocity at a given time and local coordinates. Local.
virtual scalar celerity() const
The wave celerity [m/s].
A class for handling words, derived from string.
tmp< fvMatrix< Type > > S(const Pair< tmp< volScalarField::Internal >> &, const VolField< Type > &)
defineTypeNameAndDebug(Airy, 0)
addToRunTimeSelectionTable(waveModel, Airy, dictionary)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensionedScalar cosh(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)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar cos(const dimensionedScalar &ds)