solitary.H
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-2021 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 Class
25  Foam::waveModels::solitary
26 
27 Description
28  Solitary wave model.
29 
30  Reference:
31  \verbatim
32  Dean, R. G., & Dalrymple, R. A. (1991).
33  Water wave mechanics for engineers and scientists (Vol. 2).
34  World Scientific Publishing Company.
35  \endverbatim
36 
37  See pages 314-317.
38 
39 SourceFiles
40  solitary.C
41 
42 \*---------------------------------------------------------------------------*/
43 
44 #ifndef solitary_H
45 #define solitary_H
46 
47 #include "waveModel.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 namespace waveModels
54 {
55 
56 /*---------------------------------------------------------------------------*\
57  Class solitary Declaration
58 \*---------------------------------------------------------------------------*/
59 
60 class solitary
61 :
62  public waveModel
63 {
64  // Private Data
65 
66  //- Depth [m]
67  const scalar depth_;
68 
69  //- Amplitude [m]
70  const scalar amplitude_;
71 
72  //- Offset [m]
73  const scalar offset_;
74 
75 
76  // Private Member Functions
77 
78  //- The wavenumber [1/m]
79  scalar k() const;
80 
81  //- The dimensionless amplitude [1]
82  scalar alpha() const;
83 
84  //- The wave celerity [m/s]
85  scalar celerity() const;
86 
87  //- The evolution parameter [1]
88  // This is analogous to the oscillation angle of a periodic wave
89  tmp<scalarField> parameter
90  (
91  const scalar t,
92  const scalarField& x
93  ) const;
94 
95  //- The dimensionless elevation [1]
97  (
98  const scalar t,
99  const scalarField& x
100  ) const;
101 
102 
103 public:
104 
105  //- Runtime type information
106  TypeName("solitary");
107 
108 
109  // Constructors
110 
111  //- Construct a copy
112  solitary(const solitary& wave);
113 
114  //- Construct from a dictionary and gravity
115  solitary(const dictionary& dict, const scalar g);
116 
117  //- Construct a clone
118  virtual autoPtr<waveModel> clone() const
119  {
120  return autoPtr<waveModel>(new solitary(*this));
121  }
122 
123 
124  //- Destructor
125  virtual ~solitary();
126 
127 
128  // Member Functions
129 
130  // Access
131 
132  //- Get the depth
133  scalar depth() const
134  {
135  return depth_;
136  }
137 
138  //- Get the amplitude
139  scalar amplitude() const
140  {
141  return amplitude_;
142  }
143 
144  //- Get the offset
145  scalar offset() const
146  {
147  return offset_;
148  }
149 
150  //- Get the wave elevation at a given time and local coordinates. Local
151  // x is aligned with the direction of propagation.
153  (
154  const scalar t,
155  const scalarField& x
156  ) const;
157 
158  //- Get the wave velocity at a given time and local coordinates. Local
159  // x is aligned with the direction of propagation, and z with negative
160  // gravity.
162  (
163  const scalar t,
164  const vector2DField& xz
165  ) const;
166 
167  //- Write
168  virtual void write(Ostream& os) const;
169 };
170 
171 
172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173 
174 } // End namespace waveModels
175 } // End namespace Foam
176 
177 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
178 
179 #endif
180 
181 // ************************************************************************* //
solitary(const solitary &wave)
Construct a copy.
Definition: solitary.C:86
dictionary dict
scalar depth() const
Get the depth.
Definition: solitary.H:132
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
scalar offset() const
Get the offset.
Definition: solitary.H:144
scalar g() const
Get the value of gravity.
Definition: waveModel.H:119
virtual void write(Ostream &os) const
Write.
Definition: solitary.C:154
Generic base class for waves. Derived classes must implement field functions which return the elevati...
Definition: waveModel.H:55
virtual tmp< vector2DField > velocity(const scalar t, const vector2DField &xz) const
Get the wave velocity at a given time and local coordinates. Local.
Definition: solitary.C:127
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
A class for managing temporary objects.
Definition: PtrList.H:53
virtual autoPtr< waveModel > clone() const
Construct a clone.
Definition: solitary.H:117
scalar amplitude() const
Get the amplitude.
Definition: solitary.H:138
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.
virtual ~solitary()
Destructor.
Definition: solitary.C:110
TypeName("solitary")
Runtime type information.
Namespace for OpenFOAM.
Solitary wave model.
Definition: solitary.H:59
virtual tmp< scalarField > elevation(const scalar t, const scalarField &x) const
Get the wave elevation at a given time and local coordinates. Local.
Definition: solitary.C:117