StokesIIWaveModel.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 -------------------------------------------------------------------------------
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 "StokesIIWaveModel.H"
30 #include "mathematicalConstants.H"
32 
33 using namespace Foam::constant;
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace waveModels
40 {
41  defineTypeNameAndDebug(StokesII, 0);
43  (
44  waveModel,
45  StokesII,
46  patch
47  );
48 }
49 }
50 
51 
52 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
53 
54 Foam::scalar Foam::waveModels::StokesII::eta
55 (
56  const scalar H,
57  const scalar h,
58  const scalar Kx,
59  const scalar x,
60  const scalar Ky,
61  const scalar y,
62  const scalar omega,
63  const scalar t,
64  const scalar phase
65 ) const
66 {
67  const scalar k = sqrt(Kx*Kx + Ky*Ky);
68  const scalar sigma = tanh(k*h);
69  const scalar phaseTot = Kx*x + Ky*y - omega*t + phase;
70 
71  return
72  H*0.5*cos(phaseTot)
73  + k*H*H/4.0*(3.0 - sigma*sigma)/(4.0*pow3(sigma))*cos(2.0*phaseTot);
74 }
75 
76 
77 // * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
78 
80 (
81  const scalar H,
82  const scalar h,
83  const scalar Kx,
84  const scalar x,
85  const scalar Ky,
86  const scalar y,
87  const scalar omega,
88  const scalar t,
89  const scalar phase,
90  const scalar z
91 ) const
92 {
93  const scalar k = sqrt(Kx*Kx + Ky*Ky);
94  const scalar phaseTot = Kx*x + Ky*y - omega*t + phase;
95 
96  scalar u =
97  H*0.5*omega*cos(phaseTot)*cosh(k*z)/sinh(k*h)
98  + 3.0/4.0*H*H/4.0*omega*k*cosh(2.0*k*z)/pow4(sinh(k*h))*cos(2.0*phaseTot);
99  scalar w =
100  H*0.5*omega*sin(phaseTot)*sinh(k*z)/sinh(k*h)
101  + 3.0/4.0*H*H/4.0*omega*k*sinh(2.0*k*z)/pow4(sinh(k*h))*sin(2.0*phaseTot);
102  scalar v = u*sin(waveAngle_);
103  u *= cos(waveAngle_);
104 
105  return vector(u, v, w);
106 }
107 
108 
110 (
111  const scalar t,
112  const scalar tCoeff,
113  scalarField& level
114 ) const
115 {
116  const scalar waveOmega = mathematical::twoPi/wavePeriod_;
117  const scalar waveK = mathematical::twoPi/waveLength_;
118  const scalar waveKx = waveK*cos(waveAngle_);
119  const scalar waveKy = waveK*sin(waveAngle_);
120 
121  forAll(level, paddlei)
122  {
123  const scalar eta =
124  this->eta
125  (
126  waveHeight_,
127  waterDepthRef_,
128  waveKx,
129  xPaddle_[paddlei],
130  waveKy,
131  yPaddle_[paddlei],
132  waveOmega,
133  t,
134  wavePhase_
135  );
136 
137  level[paddlei] = waterDepthRef_ + tCoeff*eta;
138  }
139 }
140 
141 
142 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
143 
145 (
146  const dictionary& dict,
147  const fvMesh& mesh,
148  const polyPatch& patch,
149  const bool readFields
150 )
151 :
152  StokesI(dict, mesh, patch, false)
153 {
154  if (readFields)
155  {
156  readDict(dict);
157  }
158 }
159 
160 
161 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
162 
164 {
165  if (StokesI::readDict(overrideDict))
166  {
167  return true;
168  }
169 
170  return false;
171 }
172 
173 
175 {
176  StokesI::info(os);
177 }
178 
179 
180 // ************************************************************************* //
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
mathematicalConstants.H
Foam::waveModels::StokesII::StokesII
StokesII(const dictionary &dict, const fvMesh &mesh, const polyPatch &patch, const bool readFields=true)
Constructor.
Definition: StokesIIWaveModel.C:145
Foam::cosh
dimensionedScalar cosh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:271
StokesIIWaveModel.H
Foam::phase
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:54
Foam::constant
Different types of constants.
Definition: atomicConstants.C:38
Foam::waveModels::StokesII::setLevel
virtual void setLevel(const scalar t, const scalar tCoeff, scalarField &level) const
Set the water level.
Definition: StokesIIWaveModel.C:110
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))
Foam::IOobject::info
InfoProxy< IOobject > info() const
Return info proxy.
Definition: IOobject.H:689
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::pow4
dimensionedScalar pow4(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:100
Foam::Field< scalar >
Foam::tanh
dimensionedScalar tanh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:272
Foam::pow3
dimensionedScalar pow3(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:89
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::StokesII::UfBase
virtual vector UfBase(const scalar H, const scalar h, const scalar Kx, const scalar x, const scalar Ky, const scalar y, const scalar omega, const scalar t, const scalar phase, const scalar z) const
Wave velocity.
Definition: StokesIIWaveModel.C:80
Foam::constant::mathematical::twoPi
constexpr scalar twoPi(2 *M_PI)
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::waveModels::StokesII::readDict
virtual bool readDict(const dictionary &overrideDict)
Read from dictionary.
Definition: StokesIIWaveModel.C:163
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::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
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
k
label k
Boltzmann constant.
Definition: LISASMDCalcMethod2.H:41
x
x
Definition: LISASMDCalcMethod2.H:52
Foam::waveModels::StokesI
Stokes I wave model.
Definition: StokesIWaveModel.H:50
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
sigma
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
y
scalar y
Definition: LISASMDCalcMethod1.H:14
Foam::cos
dimensionedScalar cos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:265
Foam::sinh
dimensionedScalar sinh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:270