angularOscillatingVelocityPointPatchVectorField.C
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) 2011-2019 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 \*---------------------------------------------------------------------------*/
25 
27 #include "pointPatchFields.H"
29 #include "Time.H"
30 #include "polyMesh.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
41 (
42  const pointPatch& p,
44 )
45 :
47  axis_(Zero),
48  origin_(Zero),
49  angle0_(0.0),
50  amplitude_(0.0),
51  omega_(0.0),
52  p0_(p.localPoints())
53 {}
54 
55 
58 (
59  const pointPatch& p,
61  const dictionary& dict
62 )
63 :
65  axis_(dict.lookup("axis")),
66  origin_(dict.lookup("origin")),
67  angle0_(dict.lookup<scalar>("angle0")),
68  amplitude_(dict.lookup<scalar>("amplitude")),
69  omega_(dict.lookup<scalar>("omega"))
70 {
71  if (!dict.found("value"))
72  {
73  updateCoeffs();
74  }
75 
76  if (dict.found("p0"))
77  {
78  p0_ = vectorField("p0", dict , p.size());
79  }
80  else
81  {
82  p0_ = p.localPoints();
83  }
84 }
85 
86 
89 (
91  const pointPatch& p,
93  const pointPatchFieldMapper& mapper
94 )
95 :
96  fixedValuePointPatchField<vector>(ptf, p, iF, mapper),
97  axis_(ptf.axis_),
98  origin_(ptf.origin_),
99  angle0_(ptf.angle0_),
100  amplitude_(ptf.amplitude_),
101  omega_(ptf.omega_),
102  p0_(ptf.p0_)
103 {}
104 
105 
108 (
111 )
112 :
114  axis_(ptf.axis_),
115  origin_(ptf.origin_),
116  angle0_(ptf.angle0_),
117  amplitude_(ptf.amplitude_),
118  omega_(ptf.omega_),
119  p0_(ptf.p0_)
120 {}
121 
122 
123 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
124 
126 (
127  const pointPatchFieldMapper& m
128 )
129 {
131 
132  m(p0_, p0_);
133 }
134 
135 
137 (
138  const pointPatchField<vector>& ptf,
139  const labelList& addr
140 )
141 {
143  refCast<const angularOscillatingVelocityPointPatchVectorField>(ptf);
144 
146 
147  p0_.rmap(aOVptf.p0_, addr);
148 }
149 
150 
152 {
153  if (this->updated())
154  {
155  return;
156  }
157 
158  const polyMesh& mesh = this->internalField().mesh()();
159  const Time& t = mesh.time();
160  const pointPatch& p = this->patch();
161 
162  scalar angle = angle0_ + amplitude_*sin(omega_*t.value());
163  vector axisHat = axis_/mag(axis_);
164  vectorField p0Rel(p0_ - origin_);
165 
167  (
168  (
169  p0_
170  + p0Rel*(cos(angle) - 1)
171  + (axisHat ^ p0Rel*sin(angle))
172  + (axisHat & p0Rel)*(1 - cos(angle))*axisHat
173  - p.localPoints()
174  )/t.deltaTValue()
175  );
176 
178 }
179 
180 
182 (
183  Ostream& os
184 ) const
185 {
187  writeEntry(os, "axis", axis_);
188  writeEntry(os, "origin", origin_);
189  writeEntry(os, "angle0", angle0_);
190  writeEntry(os, "amplitude", amplitude_);
191  writeEntry(os, "omega", omega_);
192  writeEntry(os, "p0", p0_);
193  writeEntry(os, "value", *this);
194 }
195 
196 
197 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
198 
200 (
203 );
204 
205 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
206 
207 } // End namespace Foam
208 
209 // ************************************************************************* //
dictionary dict
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:643
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
Foam::pointPatchFieldMapper.
label size() const
Return size.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Macros for easy insertion into run-time selection tables.
friend Ostream & operator(Ostream &, const Field< vector > &)
dynamicFvMesh & mesh
dimensionedScalar cos(const dimensionedScalar &ds)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
angularOscillatingVelocityPointPatchVectorField(const pointPatch &, const DimensionedField< vector, pointMesh > &)
Construct from patch and internal field.
const pointPatch & patch() const
Return patch.
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
const Type & value() const
Return const reference to value.
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
static const zero Zero
Definition: zero.H:97
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
dimensionedScalar sin(const dimensionedScalar &ds)
virtual void rmap(const pointPatchField< vector > &, const labelList &)
Reverse map the given pointPatchField onto this pointPatchField.
const Time & time() const
Return time.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
virtual const vectorField & localPoints() const =0
Return mesh points.
const DimensionedField< vector, pointMesh > & internalField() const
Return dimensioned internal field reference.
virtual void write(Ostream &) const
Write.
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:56
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
bool updated() const
Return true if the boundary condition has already been updated.
makePointPatchTypeField(pointPatchVectorField, solidBodyMotionDisplacementPointPatchVectorField)
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
volScalarField & p
virtual void rmap(const pointPatchField< Type > &, const labelList &)
Reverse map the given pointPatchField onto this pointPatchField.
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844