GrimshawWaveModel.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) 2017 IH-Cantabria
9  Copyright (C) 2017 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 
29 #include "GrimshawWaveModel.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace waveModels
37 {
38  defineTypeNameAndDebug(Grimshaw, 0);
40  (
41  waveModel,
42  Grimshaw,
43  patch
44  );
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
50 
52 (
53  const scalar H,
54  const scalar h
55 ) const
56 {
57  scalar eps = H/h;
58 
59  return sqrt(0.75*eps)*(1.0 - 0.625*eps + (71.0/128.0)*eps*eps);
60 }
61 
62 
64 (
65  const scalar H,
66  const scalar h,
67  const scalar x,
68  const scalar y,
69  const scalar theta,
70  const scalar t,
71  const scalar X0
72 ) const
73 {
74  const scalar eps = H/h;
75  const scalar eps2 = eps*eps;
76  const scalar eps3 = eps*eps2;
77 
78  const scalar C =
79  sqrt(mag(g_)*h)*sqrt(1.0 + eps - 0.05*eps2 - (3.0/70.0)*eps3);
80 
81  const scalar ts = 3.5*h/sqrt(H/h);
82  const scalar xa = -C*t + ts - X0 + x*cos(theta) + y*sin(theta);
83  const scalar alfa = this->alfa(H, h);
84 
85  const scalar s = (1.0)/(cosh(alfa*(xa/h)));
86  const scalar s2 = s*s;
87  const scalar q = tanh(alfa*(xa/h));
88  const scalar q2 = q*q;
89 
90  return
91  h
92  *(
93  eps*s2
94  - 0.75*eps2*s2*q2
95  + eps3*(0.625*s2*q2 - 1.2625*s2*s2*q2)
96  );
97 }
98 
99 
101 (
102  const scalar H,
103  const scalar h,
104  const scalar x,
105  const scalar y,
106  const scalar theta,
107  const scalar t,
108  const scalar X0,
109  const scalar z
110 ) const
111 {
112  const scalar eps = H/h;
113  const scalar eps2 = eps*eps;
114  const scalar eps3 = eps*eps2;
115 
116  const scalar C =
117  sqrt(mag(g_)*h)*sqrt(1.0 + eps - 0.05*eps2 - (3.0/70.0)*eps3);
118 
119  const scalar ts = 3.5*h/sqrt(eps);
120  const scalar xa = -C*t + ts - X0 + x*cos(theta) + y*sin(theta);
121  const scalar alfa = this->alfa(H, h);
122 
123  const scalar s = (1.0)/(cosh(alfa*(xa/h)));
124  const scalar s2 = s*s;
125  const scalar s4 = s2*s2;
126  const scalar s6 = s2*s4;
127 
128  const scalar zbyh = z/h;
129  const scalar zbyh2 = zbyh*zbyh;
130  const scalar zbyh4 = zbyh2*zbyh2;
131 
132  scalar outa = eps*s2 - eps2*(-0.25*s2 + s4 + zbyh2*(1.5*s2 - 2.25*s4));
133  scalar outb = 0.475*s2 + 0.2*s4 - 1.2*s6;
134  scalar outc = zbyh2*(-1.5*s2 - 3.75*s4 + 7.5*s6);
135  scalar outd = zbyh4*(-0.375*s2 + (45.0/16.0)*s4 - (45.0/16.0)*s6);
136 
137  scalar u = sqrt(mag(g_)*h)*(outa - eps3*(outb + outc + outd));
138 
139  outa = eps*s2 - eps2*(0.375*s2 + 2*s4 + zbyh2*(0.5*s2 - 1.5*s4));
140  outb = (49.0/640.0)*s2 - 0.85*s4 - 3.6*s6;
141  outc = zbyh2*((-13.0/16.0)*s2 -(25.0/16.0)*s4 + 7.5*s6);
142  outd = zbyh4*(-0.075*s2 -1.125*s4 - (27.0/16.0)*s6);
143 
144  const scalar w = sqrt(mag(g_)*h)*(outa - eps3*(outb + outc + outd));
145 
146  const scalar v = u*sin(waveAngle_);
147  u *= cos(waveAngle_);
148 
149  return vector(u, v, w);
150 }
151 
152 
154 (
155  const scalar t,
156  const scalar tCoeff,
157  scalarField& level
158 ) const
159 {
160  forAll(level, paddlei)
161  {
162  const scalar eta =
163  this->eta
164  (
165  waveHeight_,
166  waterDepthRef_,
167  xPaddle_[paddlei],
168  yPaddle_[paddlei],
169  waveAngle_,
170  t,
171  x0_
172  );
173 
174  level[paddlei] = waterDepthRef_ + tCoeff*eta;
175  }
176 }
177 
178 
179 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
180 
182 (
183  const dictionary& dict,
184  const fvMesh& mesh,
185  const polyPatch& patch,
186  const bool readFields
187 )
188 :
189  solitaryWaveModel(dict, mesh, patch, false)
190 {
191  if (readFields)
192  {
193  readDict(dict);
194  }
195 }
196 
197 
198 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
199 
201 {
202  if (solitaryWaveModel::readDict(overrideDict))
203  {
204  return true;
205  }
206 
207  return false;
208 }
209 
210 
212 (
213  const scalar t,
214  const scalar tCoeff,
215  const scalarField& level
216 )
217 {
218  forAll(U_, facei)
219  {
220  // Fraction of geometry represented by paddle - to be set
221  scalar fraction = 1;
222 
223  // Height - to be set
224  scalar z = 0;
225 
226  setPaddlePropeties(level, facei, fraction, z);
227 
228  if (fraction > 0)
229  {
230  const label paddlei = faceToPaddle_[facei];
231 
232  const vector Uf = this->Uf
233  (
234  waveHeight_,
235  waterDepthRef_,
236  xPaddle_[paddlei],
237  yPaddle_[paddlei],
238  waveAngle_,
239  t,
240  x0_,
241  z
242  );
243 
244  U_[facei] = fraction*Uf*tCoeff;
245  }
246  }
247 }
248 
249 
251 {
253 }
254 
255 
256 // ************************************************************************* //
Foam::waveModels::solitaryWaveModel
Definition: solitaryWaveModel.H:49
Foam::waveModels::Grimshaw::setLevel
virtual void setLevel(const scalar t, const scalar tCoeff, scalarField &level) const
Set the water level.
Definition: GrimshawWaveModel.C:154
Foam::waveModels::Grimshaw::Uf
virtual vector Uf(const scalar H, const scalar h, const scalar x, const scalar y, const scalar theta, const scalar t, const scalar X0, const scalar z) const
Wave velocity.
Definition: GrimshawWaveModel.C:101
s
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputSpray.H:25
Foam::cosh
dimensionedScalar cosh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:271
X0
scalarList X0(nSpecie, Zero)
Foam::sin
dimensionedScalar sin(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:264
H
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
Uf
autoPtr< surfaceVectorField > Uf
Definition: createUfIfPresent.H:33
Foam::waveModels::addToRunTimeSelectionTable
addToRunTimeSelectionTable(waveModel, shallowWaterAbsorption, patch)
Foam::waveModels::Grimshaw::readDict
virtual bool readDict(const dictionary &overrideDict)
Read from dictionary.
Definition: GrimshawWaveModel.C:200
Foam::waveModels::solitaryWaveModel::readDict
virtual bool readDict(const dictionary &overrideDict)
Read from dictionary.
Definition: solitaryWaveModel.C:87
Foam::IOobject::info
InfoProxy< IOobject > info() const
Return info proxy.
Definition: IOobject.H:689
Foam::waveModels::Grimshaw::alfa
virtual scalar alfa(const scalar H, const scalar h) const
Definition: GrimshawWaveModel.C:52
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::waveModels::Grimshaw::setVelocity
virtual void setVelocity(const scalar t, const scalar tCoeff, const scalarField &level)
Calculate the wave model velocity.
Definition: GrimshawWaveModel.C:212
Foam::Field< scalar >
Foam::tanh
dimensionedScalar tanh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:272
Foam::constant::universal::h
const dimensionedScalar h
Planck constant.
Definition: setRegionSolidFields.H:33
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::waveModels::Grimshaw::eta
virtual scalar eta(const scalar H, const scalar h, const scalar x, const scalar y, const scalar theta, const scalar t, const scalar X0) const
Wave height.
Definition: GrimshawWaveModel.C:64
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
Foam::waveModels::Grimshaw::Grimshaw
Grimshaw(const dictionary &dict, const fvMesh &mesh, const polyPatch &patch, const bool readFields=true)
Constructor.
Definition: GrimshawWaveModel.C:182
Foam::readFields
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
Definition: ReadFieldsTemplates.C:312
GrimshawWaveModel.H
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::Vector< scalar >
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
x
x
Definition: LISASMDCalcMethod2.H:52
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::C
Graphite solid properties.
Definition: C.H:50
Foam::waveModels::defineTypeNameAndDebug
defineTypeNameAndDebug(waveAbsorptionModel, 0)
y
scalar y
Definition: LISASMDCalcMethod1.H:14
Foam::cos
dimensionedScalar cos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:265