46 const scalar amplitude,
51 static const scalar kdGreat =
log(great);
77 const word& modelName,
78 scalar (*modelCelerity)(scalar, scalar, scalar, scalar)
99 static const scalar kdGreat =
log(great);
100 const scalar kd =
min(
max(
k()*depth(), - kdGreat), kdGreat);
101 const scalar ka =
k()*amplitude(t);
103 const scalar
S = deep() ? 0 : 1/
cosh(2*kd),
T = deep() ? 1 :
tanh(kd);
124 9.0/128/(3 + 2*
S)/(4 +
S)/
pow6(1 -
S)
131 5.0/384/(3 + 2*
S)/(4 +
S)/
pow6(1 -
S)
140 <<
"B42 = " << B42 <<
endl
141 <<
"B44 = " << B44 <<
endl
142 <<
"B53 = " << B53 <<
endl
143 <<
"B55 = " << B55 <<
endl;
154 +
pow5(ka)*(- (B53 + B55)*
cos(phi) + B53*
cos(3*phi) + B55*
cos(5*phi))
165 static const scalar kdGreat =
log(great);
166 const scalar kd =
min(
max(
k()*depth(), - kdGreat), kdGreat);
167 const scalar ka =
k()*amplitude(t);
169 const scalar
S = deep() ? 0 : 1/
cosh(2*kd);
170 const scalar SByA11 = deep() ? 0 :
S*
sinh(kd);
172 const scalar A31ByA11 =
178 const scalar A33ByA11 =
184 const scalar A42ByA11 =
185 SByA11/24/
pow5(1 -
S)
190 const scalar A44ByA11 =
191 SByA11/48/(3 + 2*
S)/
pow5(1 -
S)
196 const scalar A51ByA11 =
197 1.0/64/(3 + 2*
S)/(4 +
S)/
pow6(1 -
S)
203 const scalar A53ByA11 =
204 1.0/32/(3 + 2*
S)/
pow6(1 -
S)
210 const scalar A55ByA11 =
211 1.0/64/(3 + 2*
S)/(4 +
S)/
pow6(1 -
S)
219 const scalar A11 = 1/
sinh(kd);
220 Info<<
"A31 = " << A31ByA11*A11 <<
endl
221 <<
"A33 = " << A33ByA11*A11 <<
endl
222 <<
"A42 = " << A42ByA11*A11 <<
endl
223 <<
"A44 = " << A44ByA11*A11 <<
endl
224 <<
"A51 = " << A51ByA11*A11 <<
endl
225 <<
"A53 = " << A53ByA11*A11 <<
endl
226 <<
"A55 = " << A55ByA11*A11 <<
endl;
235 pow3(ka)*(A31ByA11*v1 + A33ByA11*v3)
236 +
pow4(ka)*(A42ByA11*vi(2, t, xz) + A44ByA11*vi(4, t, xz))
237 +
pow5(ka)*(A51ByA11*v1 + A53ByA11*v3 + A55ByA11*vi(5, t, xz))
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.
scalar k() const
The angular wavenumber [rad/m].
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 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].
virtual ~Stokes5()
Destructor.
Stokes5(const dictionary &dict, const scalar g, const word &modelName=Stokes5::typeName, scalar(*modelCelerity)(scalar, scalar, scalar, scalar)=&Stokes5::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)
dimensionedScalar pow6(const dimensionedScalar &ds)
dimensionedScalar pow5(const dimensionedScalar &ds)
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)
dimensionedScalar pow4(const dimensionedScalar &ds)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar cos(const dimensionedScalar &ds)