irregularMultiDirectionalWaveModel.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015 IH-Cantabria
9 Copyright (C) 2016-2017 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
30#include "unitConversion.H"
32
33using namespace Foam::constant;
34
35// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36
37namespace Foam
38{
39namespace waveModels
40{
43 (
46 patch
47 );
48}
49}
50
51
52// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
53
55(
56 const scalar H,
57 const scalar Kx,
58 const scalar x,
59 const scalar Ky,
60 const scalar y,
61 const scalar omega,
62 const scalar t,
63 const scalar phase
64) const
65{
66 scalar phaseTot = Kx*x + Ky*y - omega*t + phase;
67 return H*0.5*cos(phaseTot);
68}
69
70
71// * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
72
74(
75 const scalar h,
76 const scalar T
77) const
78{
79 scalar L0 = mag(g_)*T*T/(2.0*mathematical::pi);
80 scalar L = L0;
81
82 for (label i=1; i<=100; ++i)
83 {
84 L = L0*tanh(2.0*mathematical::pi*h/L);
85 }
86
87 return L;
88}
89
90
92(
93 const scalar h,
94 const scalar x,
95 const scalar y,
96 const scalar t,
97 const scalar z
98) const
99{
100 scalar u = 0.0;
101 scalar v = 0.0;
102 scalar w = 0.0;
103
104 forAll(irregWaveHeights_, ii)
105 {
106 forAll(irregWaveHeights_[ii], jj)
107 {
108 scalar waveKs = mathematical::twoPi/irregWaveLengths_[ii][jj];
109 scalar waveOmegas = mathematical::twoPi/irregWavePeriods_[ii][jj];
110
111 scalar phaseTot =
112 waveKs*x*cos(irregWaveDirs_[ii][jj])
113 + waveKs*y*sin(irregWaveDirs_[ii][jj])
114 - waveOmegas*t
115 + irregWavePhases_[ii][jj];
116
117 const vector Uf = uMultiDirec
118 (
119 irregWaveHeights_[ii][jj],
120 waveOmegas,
121 phaseTot,
122 waveKs,
123 z,
124 h,
125 irregWaveDirs_[ii][jj]
126 );
127
128 u += Uf[0];
129 v += Uf[1];
130 w += Uf[2];
131 }
132 }
133
134 return vector(u, v, w);
135}
136
137
139(
140 const scalar t,
141 const scalar tCoeff,
142 scalarField& level
143) const
144{
145 forAll(level, paddlei)
146 {
147 scalar eta = 0;
148
149 forAll(irregWaveHeights_, ii)
150 {
151 forAll(irregWaveHeights_[ii], jj)
152 {
153 scalar waveKs = mathematical::twoPi/irregWaveLengths_[ii][jj];
154 scalar waveOmegas =
155 mathematical::twoPi/irregWavePeriods_[ii][jj];
156
157 eta +=
158 this->eta
159 (
160 irregWaveHeights_[ii][jj],
161 waveKs*cos(irregWaveDirs_[ii][jj]),
162 xPaddle_[paddlei],
163 waveKs*sin(irregWaveDirs_[ii][jj]),
164 yPaddle_[paddlei],
165 waveOmegas,
166 t,
167 irregWavePhases_[ii][jj]
168 );
169 }
170 }
171
172 level[paddlei] = waterDepthRef_ + tCoeff*eta;
173 }
174}
175
176
178(
179 const scalar irregH,
180 const scalar irregWaveOmega,
181 const scalar pha,
182 const scalar irregWaveKs,
183 const scalar zz,
184 const scalar hh,
185 const scalar irregDir
186) const
187{
188 const scalar ksh = irregWaveKs*hh;
189 const scalar ksz = irregWaveKs*zz;
190
191 scalar u =
192 irregH*0.5*irregWaveOmega*cos(pha)*(cosh(ksz)/sinh(ksh))*cos(irregDir);
193
194 scalar v =
195 irregH*0.5*irregWaveOmega*cos(pha)*(cosh(ksz)/sinh(ksh))*sin(irregDir);
196
197 scalar w =
198 irregH*0.5*irregWaveOmega*sin(pha)*(sinh(ksz)/sinh(ksh));
199
200 return vector(u, v, w);
201}
202
203
205(
206 const scalar t,
207 const scalar tCoeff,
208 const scalarField& level
209)
210{
211 forAll(U_, facei)
212 {
213 // Fraction of geometry represented by paddle - to be set
214 scalar fraction = 1;
215
216 // Height - to be set
217 scalar z = 0;
218
219 setPaddlePropeties(level, facei, fraction, z);
220
221 if (fraction > 0)
222 {
223 const label paddlei = faceToPaddle_[facei];
224
225 const vector Uf = this->Uf
226 (
227 waterDepthRef_,
228 xPaddle_[paddlei],
229 yPaddle_[paddlei],
230 t,
231 z
232 );
233
234 U_[facei] = fraction*Uf*tCoeff;
235 }
236 }
237}
238
239
240// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
241
243(
244 const dictionary& dict,
245 const fvMesh& mesh,
246 const polyPatch& patch,
247 const bool readFields
248)
249:
250 irregularWaveModel(dict, mesh, patch, false)
251{
252 if (readFields)
253 {
254 readDict(dict);
255 }
256}
257
258
259// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
260
262(
263 const dictionary& overrideDict
264)
265{
266 if (irregularWaveModel::readDict(overrideDict))
267 {
268 readEntry("wavePeriods", irregWavePeriods_);
269 readEntry("waveHeights", irregWaveHeights_);
270 readEntry("wavePhases", irregWavePhases_);
271 readEntry("waveDirs", irregWaveDirs_);
272
273 irregWaveLengths_ = irregWaveHeights_;
274
275 forAll(irregWaveHeights_, ii)
276 {
277 forAll(irregWaveHeights_[ii], jj)
278 {
279 irregWaveLengths_[ii][jj] =
280 waveLength(waterDepthRef_, irregWavePeriods_[ii][jj]);
281 irregWaveDirs_[ii][jj] =
282 degToRad(irregWaveDirs_[ii][jj]);
283 }
284 }
285
286 return true;
287 }
288
289 return false;
290}
291
292
294{
296
297 os << " Wave periods : " << irregWavePeriods_.size() << nl
298 << " Wave heights : " << irregWaveHeights_.size() << nl
299 << " Wave phases : " << irregWavePhases_.size() << nl
300 << " Wave lengths : " << irregWaveLengths_.size() << nl
301 << " Wave directions : " << irregWaveDirs_.size() << nl;
302}
303
304
305// ************************************************************************* //
scalar y
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
InfoProxy< ensightCells > info() const
Return info proxy.
Definition: ensightCells.H:254
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:57
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
Base class for waveModels.
Definition: waveModel.H:61
virtual void setVelocity(const scalar t, const scalar tCoeff, const scalarField &level)
Calculate the wave model velocity.
virtual vector uMultiDirec(const scalar irregH, const scalar irregWaveOmega, const scalar phaseTot, const scalar irregWaveKs, const scalar z, const scalar h, const scalar irregDir) const
virtual void setLevel(const scalar t, const scalar tCoeff, scalarField &level) const
Set the water level.
virtual scalar waveLength(const scalar h, const scalar T) const
Return the wavelength.
virtual bool readDict(const dictionary &overrideDict)
Read from dictionary.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
const volScalarField & T
dynamicFvMesh & mesh
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
autoPtr< surfaceVectorField > Uf
OBJstream os(runTime.globalPath()/outputName)
constexpr scalar pi(M_PI)
constexpr scalar twoPi(2 *M_PI)
Different types of constants.
Namespace for OpenFOAM.
dimensionedScalar cosh(const dimensionedScalar &ds)
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar tanh(const dimensionedScalar &ds)
dimensionedScalar sinh(const dimensionedScalar &ds)
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
constexpr scalar degToRad() noexcept
Multiplication factor for degrees to radians conversion.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Vector< scalar > vector
Definition: vector.H:61
dimensionedScalar cos(const dimensionedScalar &ds)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
volScalarField & h
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Unit conversion functions.
const vector L(dict.get< vector >("L"))