angularOscillatingDisplacementPointPatchVectorField.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) 2011-2016 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
29#include "pointPatchFields.H"
31#include "Time.H"
32#include "polyMesh.H"
33
34// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35
36namespace Foam
37{
38
39// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40
41angularOscillatingDisplacementPointPatchVectorField::
42angularOscillatingDisplacementPointPatchVectorField
43(
44 const pointPatch& p,
46)
47:
49 axis_(Zero),
50 origin_(Zero),
51 angle0_(0.0),
52 amplitude_(0.0),
53 omega_(0.0),
54 p0_(p.localPoints())
55{}
56
57
58angularOscillatingDisplacementPointPatchVectorField::
59angularOscillatingDisplacementPointPatchVectorField
60(
61 const pointPatch& p,
63 const dictionary& dict
64)
65:
67 axis_(dict.lookup("axis")),
68 origin_(dict.lookup("origin")),
69 angle0_(dict.get<scalar>("angle0")),
70 amplitude_(dict.get<scalar>("amplitude")),
71 omega_(dict.get<scalar>("omega"))
72{
73 if (!dict.found("value"))
74 {
76 }
77
78 if (dict.found("p0"))
79 {
80 p0_ = vectorField("p0", dict , p.size());
81 }
82 else
83 {
84 p0_ = p.localPoints();
85 }
86}
87
88
89angularOscillatingDisplacementPointPatchVectorField::
90angularOscillatingDisplacementPointPatchVectorField
91(
93 const pointPatch& p,
95 const pointPatchFieldMapper& mapper
96)
97:
98 fixedValuePointPatchField<vector>(ptf, p, iF, mapper),
99 axis_(ptf.axis_),
100 origin_(ptf.origin_),
101 angle0_(ptf.angle0_),
102 amplitude_(ptf.amplitude_),
103 omega_(ptf.omega_),
104 p0_(ptf.p0_, mapper)
105{}
106
107
108angularOscillatingDisplacementPointPatchVectorField::
109angularOscillatingDisplacementPointPatchVectorField
110(
113)
114:
116 axis_(ptf.axis_),
117 origin_(ptf.origin_),
118 angle0_(ptf.angle0_),
119 amplitude_(ptf.amplitude_),
120 omega_(ptf.omega_),
121 p0_(ptf.p0_)
122{}
123
124
125// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
126
128(
129 const pointPatchFieldMapper& m
130)
131{
133
134 p0_.autoMap(m);
135}
136
137
139(
140 const pointPatchField<vector>& ptf,
141 const labelList& addr
142)
143{
145 refCast<const angularOscillatingDisplacementPointPatchVectorField>(ptf);
146
148
149 p0_.rmap(aODptf.p0_, addr);
150}
151
152
154{
155 if (this->updated())
156 {
157 return;
158 }
159
160 const polyMesh& mesh = this->internalField().mesh()();
161 const Time& t = mesh.time();
162
163 scalar angle = angle0_ + amplitude_*sin(omega_*t.value());
164 vector axisHat = axis_/mag(axis_);
165 vectorField p0Rel(p0_ - origin_);
166
168 (
169 p0Rel*(cos(angle) - 1)
170 + (axisHat ^ p0Rel*sin(angle))
171 + (axisHat & p0Rel)*(1 - cos(angle))*axisHat
172 );
173
175}
176
177
179(
180 Ostream& os
181) const
182{
184 os.writeEntry("axis", axis_);
185 os.writeEntry("origin", origin_);
186 os.writeEntry("angle0", angle0_);
187 os.writeEntry("amplitude", amplitude_);
188 os.writeEntry("omega", omega_);
189 p0_.writeEntry("p0", os);
190 writeEntry("value", os);
191}
192
193
194// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
195
197(
200);
201
202// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
203
204} // End namespace Foam
205
206// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Mesh & mesh() const
Return mesh.
friend Ostream & operator(Ostream &, const Field< vector > &)
void autoMap(const FieldMapper &map, const bool applyFlip=true)
Map from self.
Definition: Field.C:403
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
Definition: Field.C:614
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:466
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:239
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void rmap(const pointPatchField< vector > &, const labelList &)
Reverse map the given pointPatchField onto this pointPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
void rmap(const atmBoundaryLayer &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Definition: dictionaryI.H:87
const Type & value() const
Return const reference to value.
A FixedValue boundary condition for pointField.
virtual bool write()
Write the output fields.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
Foam::pointPatchFieldMapper.
Abstract base class for point-mesh patch fields.
const DimensionedField< Type, pointMesh > & internalField() const
Return dimensioned internal field reference.
bool updated() const
Return true if the boundary condition has already been updated.
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:64
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
Lookup type of boundary radiation properties.
Definition: lookup.H:66
volScalarField & p
dynamicFvMesh & mesh
OBJstream os(runTime.globalPath()/outputName)
Namespace for OpenFOAM.
dimensionedScalar sin(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Field< vector > vectorField
Specialisation of Field<T> for vector.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
dimensionedScalar cos(const dimensionedScalar &ds)
#define makePointPatchTypeField(PatchTypeField, typePatchTypeField)
dictionary dict