Go to the documentation of this file.
45 { motionTypes::piston,
"piston" },
46 { motionTypes::flap,
"flap" },
47 { motionTypes::solitary,
"solitary" }
60 <<
"Gravity vector is not set. Please update "
61 << gf.uniformDimensionedVectorField::path()
78 for (label i=1; i<=100; ++i)
92 return max(0,
min(t/rampTime_, 1));
105 const scalar yMin = bb.
min().
y();
106 const scalar yMax = bb.
max().
y();
107 zSpan_ = bb.
max().
z() - bb.
min().
z();
109 zMinGb_ = bb.
min().
z();
113 xPaddle_.setSize(nPaddle_, 0);
114 yPaddle_.setSize(nPaddle_, 0);
116 const scalar paddleDy = (yMax - yMin)/scalar(nPaddle_);
118 for (label paddlei = 0; paddlei < nPaddle_; ++paddlei)
120 xPaddle_[paddlei] = xMid;
121 yPaddle_[paddlei] = paddlei*paddleDy + yMin + 0.5*paddleDy;
125 x_ = this->
patch().localPoints().component(0);
126 y_ = this->
patch().localPoints().component(1);
127 z_ = this->
patch().localPoints().component(2);
130 pointToPaddle_.setSize(this->
patch().size(), -1);
132 forAll(pointToPaddle_, ppi)
134 pointToPaddle_[ppi] = floor((y_[ppi] - yMin)/(paddleDy+0.01*paddleDy));
147 motionType_(motionTypes::piston),
170 motionType_(motionTypeNames.get(
"motionType",
dict)),
173 initialDepth_(
dict.
get<scalar>(
"initialDepth")),
174 wavePeriod_(
dict.
get<scalar>(
"wavePeriod")),
175 waveHeight_(
dict.
get<scalar>(
"waveHeight")),
176 wavePhase_(
dict.
get<scalar>(
"wavePhase")),
186 rampTime_(
dict.
get<scalar>(
"rampTime")),
194 <<
"Patch normal direction vector is not set. 'n' = " << n_
199 gHat_ = (
g() - n_*(n_&
g()));
200 if (
mag(gHat_) < SMALL)
203 <<
"Patch normal and gravity directions must not be aligned. "
204 <<
"'n' = " << n_ <<
" 'g' = " <<
g()
211 initialiseGeometry();
213 waterDepthRef_.setSize(nPaddle_, -1);
280 if (initialDepth_ != 0 )
282 forAll(waterDepthRef_, paddlei)
284 waterDepthRef_[paddlei] = initialDepth_;
290 <<
"initialDepth is not set. Please update "
295 Info<<
" WaterDepth at the wavepaddles = " << waterDepthRef_ <<
endl;
299 const scalar t = db().time().value() - startTime_;
309 waveLength_[padddlei] =
310 waveLength(waterDepthRef_[padddlei], wavePeriod_);
313 waveKx[padddlei] = waveK[padddlei]*
cos(waveAngle_);
314 waveKy[padddlei] = waveK[padddlei]*
sin(waveAngle_);
320 case motionTypes::flap:
327 const label paddlei = pointToPaddle_[pointi];
329 const scalar phaseTot =
330 waveKx[paddlei]*xPaddle_[paddlei]
331 + waveKy[paddlei]*yPaddle_[paddlei];
333 const scalar depthRef = waterDepthRef_[paddlei];
334 const scalar kh = waveK[paddlei]*depthRef;
341 const scalar boardStroke = waveHeight_/m1;
343 motionX[pointi] = 0.5*boardStroke*
sin(phaseTot -
sigma*t);
348 sqr(waveHeight_)/(16*depthRef)
354 motionX[pointi] *= 1.0 + (pz - zMinGb_ - depthRef)/depthRef;
362 case motionTypes::piston:
369 const label paddlei = pointToPaddle_[pointi];
371 const scalar phaseTot =
372 waveKx[paddlei]*xPaddle_[paddlei]
373 + waveKy[paddlei]*yPaddle_[paddlei];
375 const scalar depthRef = waterDepthRef_[paddlei];
376 const scalar kh = waveK[paddlei]*depthRef;
377 const scalar m1 = 2*(
cosh(2*kh) - 1.0)/(
sinh(2*kh) + 2*kh);
378 const scalar boardStroke = waveHeight_/m1;
380 motionX[pointi] = 0.5*boardStroke*
sin(phaseTot -
sigma*t);
395 case motionTypes::solitary:
399 const scalar magG =
mag(
g());
403 const label paddlei = pointToPaddle_[pointi];
404 const scalar depthRef = waterDepthRef_[paddlei];
407 const scalar celerity =
sqrt(magG*(depthRef + waveHeight_));
408 const scalar stroke =
sqrt(16*waveHeight_*depthRef/3.0);
409 const scalar
hr = waveHeight_/depthRef;
410 wavePeriod_ = 2.0/(
kappa*celerity)*(3.8 +
hr);
411 const scalar tSolitary = -0.5*wavePeriod_ + t;
417 const scalar
error = 0.001;
422 - (theta1 -
kappa*celerity*tSolitary +
hr*
tanh(theta1))
423 /(1.0 +
hr*(1.0/
cosh(theta1))*(1.0/
cosh(theta1)));
425 er =
mag(theta1 - theta2);
430 waveHeight_/(
kappa*depthRef)*
tanh(theta1) + 0.5*stroke;
440 <<
"Unhandled enumeration " << motionTypeNames[motionType_]
452 os.
writeEntry(
"motionType", motionTypeNames[motionType_]);
463 writeEntry(
"value", os);
const vector L(dict.get< vector >("L"))
const Cmpt & x() const
Access to the vector x component.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
scalar waveAngle_
Wave angle.
label nPaddle_
Number of wave paddles.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
vector gHat_
Vertical direction.
dimensionedScalar cosh(const dimensionedScalar &ds)
scalar startTime_
Start time.
static constexpr const zero Zero
Global zero (0)
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
scalar waveHeight_
Wave height.
dimensionedScalar sin(const dimensionedScalar &ds)
const objectRegistry & db() const
Return local objectRegistry.
scalar secondOrder_
On/off second order calculation switch.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
const point & max() const
Maximum describing the bounding box.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Basic pointPatch represents a set of points from the mesh.
scalar wavePhase_
Wave phase.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Abstract base class for point-mesh patch fields.
const Cmpt & z() const
Access to the vector z component.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
volScalarField & h
Planck constant.
Foam::pointPatchFieldMapper.
const point & min() const
Minimum describing the bounding box.
scalar wavePeriod_
Wave period.
#define forAll(list, i)
Loop across all elements in list.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
const dimensionedScalar hr
Reduced Planck constant: default SI units: [J/s].
makePointPatchTypeField(pointPatchVectorField, solidBodyMotionDisplacementPointPatchVectorField)
Point motion boundary condition to generate waves based on either piston or flap motions.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
dimensionedScalar tanh(const dimensionedScalar &ds)
static const Enum< motionTypes > motionTypeNames
Names for motion types.
dimensionedScalar pow3(const dimensionedScalar &ds)
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual void write(Ostream &) const
Write.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
waveMakerPointPatchVectorField(const pointPatch &, const DimensionedField< vector, pointMesh > &)
Construct from patch and internal field.
Gravitational acceleration vector Although termed a MeshObject it is registered on Time only and thus...
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
constexpr scalar twoPi(2 *M_PI)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
motionTypes motionType_
Motion type.
const uniformDimensionedVectorField & g
virtual scalar waveLength(const scalar h, const scalar T)
Dispersion equation.
errorManip< error > abort(error &err)
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
void operator=(const Field< Type > &)
Copy assignment.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const Cmpt & y() const
Access to the vector y component.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr scalar pi(M_PI)
const std::string patch
OpenFOAM patch number as a std::string.
scalar rampTime_
Ramp time.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
virtual void write(Ostream &) const
Write.
const volScalarField & Cp
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
A bounding box defined in terms of min/max extrema points.
vector n_
Patch normal direction.
virtual void initialiseGeometry()
Initialise.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
scalar initialDepth_
Initial water depth.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
virtual scalar timeCoeff(const scalar t) const
Return the time scaling coefficient.
Class to handle errors and exceptions in a simple, consistent stream-based manner.
const vector & g()
Return the gravitational acceleration.
static const gravity & New(const Time &runTime)
Construct on Time.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar sinh(const dimensionedScalar &ds)