Airy.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "Airy.H"
27 #include "mathematicalConstants.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace waveModels
35 {
38 }
39 }
40 
41 
42 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
43 
44 Foam::scalar Foam::waveModels::Airy::k() const
45 {
46  return 2*Foam::constant::mathematical::pi/length_;
47 }
48 
49 
51 {
52  return sqrt(g()/k()*tanh(k()*depth()));
53 }
54 
55 
57 (
58  const scalar t,
59  const scalarField& x
60 ) const
61 {
62  return phase_ + k()*(x - celerity()*t);
63 }
64 
65 
67 {
68  return k()*depth() > log(great);
69 }
70 
71 
73 (
74  const label i,
75  const scalar t,
76  const vector2DField& xz
77 ) const
78 {
79  const scalarField x(xz.component(0));
80  const scalarField z(xz.component(1));
81 
82  const scalarField phi(angle(t, x));
83 
84  const scalar kzGreat = log(i*great);
85  const scalarField kz(min(max(k()*z, - kzGreat), kzGreat));
86 
87  if (deep())
88  {
89  return i*exp(kz)*zip(cos(i*phi), sin(i*phi));
90  }
91  else
92  {
93  const scalar kd = k()*depth();
94  const scalarField kdz(max(scalar(0), kd + kz));
95  return i*zip(cosh(i*kdz)*cos(i*phi), sinh(i*kdz)*sin(i*phi))/sinh(kd);
96  }
97 }
98 
99 
100 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
101 
103 :
104  waveModel(wave),
105  length_(wave.length_),
106  phase_(wave.phase_),
107  depth_(wave.depth_)
108 {}
109 
110 
112 (
113  const objectRegistry& db,
114  const dictionary& dict
115 )
116 :
117  waveModel(db, dict),
118  length_(dict.lookup<scalar>("length")),
119  phase_(dict.lookup<scalar>("phase")),
120  depth_(dict.lookupOrDefault<scalar>("depth", log(2*great)/k()))
121 {}
122 
123 
124 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
125 
127 {}
128 
129 
130 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
131 
133 (
134  const scalar t,
135  const scalarField& x
136 ) const
137 {
138  return amplitude(t)*cos(angle(t, x));
139 }
140 
141 
143 (
144  const scalar t,
145  const vector2DField& xz
146 ) const
147 {
148  const scalar ka = k()*amplitude(t);
149 
150  return celerity()*ka*vi(1, t, xz);
151 }
152 
153 
155 (
156  const scalar t,
157  const vector2DField& xz
158 ) const
159 {
160  // It is a fluke of the formulation that the time derivative of the velocity
161  // potential equals the x-derivative multiplied by the celerity. This allows
162  // for this shortcut in evaluating the unsteady pressure.
163  return celerity()*velocity(t, xz)->component(0);
164 }
165 
166 
168 {
169  waveModel::write(os);
170 
171  writeEntry(os, "length", length_);
172  writeEntry(os, "phase", phase_);
173  if (!deep())
174  {
175  writeEntry(os, "depth", depth_);
176  }
177 }
178 
179 
180 // ************************************************************************* //
tmp< vector2DField > zip(const tmp< scalarField > &x, const tmp< scalarField > &y)
Definition: vector2DField.C:31
dimensionedScalar tanh(const dimensionedScalar &ds)
Airy(const Airy &wave)
Construct a copy.
Definition: Airy.C:102
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
dimensionedScalar log(const dimensionedScalar &ds)
scalar celerity() const
The wave celerity [m/s].
Definition: Airy.C:50
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
defineTypeNameAndDebug(Airy, 0)
virtual tmp< scalarField > elevation(const scalar t, const scalarField &x) const
Get the wave elevation at a given time and local coordinates. Local.
Definition: Airy.C:133
dimensionedScalar sqrt(const dimensionedScalar &ds)
First-order wave model.
Definition: Airy.H:59
bool deep() const
Return whether shallow and intermediate effects are to be omitted.
Definition: Airy.C:66
virtual tmp< vector2DField > velocity(const scalar t, const vector2DField &xz) const
Get the wave velocity at a given time and local coordinates. Local.
Definition: Airy.C:143
label k
Boltzmann constant.
Macros for easy insertion into run-time selection tables.
virtual tmp< scalarField > pressure(const scalar t, const vector2DField &xz) const
Get the wave pressure at a given time and local coordinates. Local.
Definition: Airy.C:155
tmp< scalarField > angle(const scalar t, const scalarField &x) const
Angle of the oscillation [rad].
Definition: Airy.C:57
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
Generic base class for waves. Derived classes must implement field functions which return the elevati...
Definition: waveModel.H:56
scalar k() const
The angular wavenumber [rad/m].
Definition: Airy.C:44
phi
Definition: correctPhi.H:3
addToRunTimeSelectionTable(waveModel, Airy, objectRegistry)
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:441
waveModel(const waveModel &wave)
Construct a copy.
Definition: waveModel.C:41
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
dimensionedScalar sin(const dimensionedScalar &ds)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
dimensionedScalar sinh(const dimensionedScalar &ds)
tmp< vector2DField > vi(const label i, const scalar t, const vector2DField &xz) const
Return the non-dimensionalised i-th harmonic of the velocity.
Definition: Airy.C:73
scalar amplitude(const scalar t) const
Get the amplitude.
Definition: waveModel.H:123
virtual void write(Ostream &os) const
Write.
Definition: Airy.C:167
dimensionedScalar cosh(const dimensionedScalar &ds)
A class for managing temporary objects.
Definition: PtrList.H:53
Registry of regIOobjects.
const dimensionedVector & g
virtual ~Airy()
Destructor.
Definition: Airy.C:126
virtual void write(Ostream &os) const
Write.
Definition: waveModel.C:73
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844