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 -------------------------------------------------------------------------------
10 License
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 
36 namespace Foam
37 {
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
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 
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  {
75  updateCoeffs();
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 
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 
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 // ************************************************************************* //
Foam::angularOscillatingDisplacementPointPatchVectorField
Foam::angularOscillatingDisplacementPointPatchVectorField.
Definition: angularOscillatingDisplacementPointPatchVectorField.H:50
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::dictionary::found
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
Foam::sin
dimensionedScalar sin(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:264
Foam::Field< vector >::operator
friend Ostream & operator(Ostream &, const Field< vector > &)
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::dictionary::get
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:107
Foam::pointPatch
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:58
polyMesh.H
Foam::angularOscillatingDisplacementPointPatchVectorField::write
virtual void write(Ostream &) const
Write.
Definition: angularOscillatingDisplacementPointPatchVectorField.C:179
Foam::pointPatchField< vector >
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::pointPatchFieldMapper
Foam::pointPatchFieldMapper.
Definition: pointPatchFieldMapper.H:48
Foam::vectorField
Field< vector > vectorField
Specialisation of Field<T> for vector.
Definition: primitiveFieldsFwd.H:54
angularOscillatingDisplacementPointPatchVectorField.H
Foam::fixedValuePointPatchField< vector >
Foam::makePointPatchTypeField
makePointPatchTypeField(pointPatchVectorField, solidBodyMotionDisplacementPointPatchVectorField)
Foam::Field< vector >
Foam::valuePointPatchField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: valuePointPatchField.C:135
Foam::pointPatchField< vector >::internalField
const DimensionedField< vector, pointMesh > & internalField() const
Return dimensioned internal field reference.
Definition: pointPatchField.H:275
Foam::pointPatchField::write
virtual void write(Ostream &) const
Write.
Definition: pointPatchField.C:117
Foam::angularOscillatingDisplacementPointPatchVectorField::autoMap
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: angularOscillatingDisplacementPointPatchVectorField.C:128
Foam::angularOscillatingDisplacementPointPatchVectorField::rmap
virtual void rmap(const pointPatchField< vector > &, const labelList &)
Reverse map the given pointPatchField onto this pointPatchField.
Definition: angularOscillatingDisplacementPointPatchVectorField.C:139
Foam::dictionary::lookup
ITstream & lookup(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionary.C:386
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
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::valuePointPatchField::autoMap
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: valuePointPatchField.C:108
Foam::pointPatchField< vector >::updated
bool updated() const
Return true if the boundary condition has already been updated.
Definition: pointPatchField.H:311
Foam::angularOscillatingDisplacementPointPatchVectorField::angularOscillatingDisplacementPointPatchVectorField
angularOscillatingDisplacementPointPatchVectorField(const pointPatch &, const DimensionedField< vector, pointMesh > &)
Construct from patch and internal field.
Definition: angularOscillatingDisplacementPointPatchVectorField.C:43
Time.H
Foam::angularOscillatingDisplacementPointPatchVectorField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: angularOscillatingDisplacementPointPatchVectorField.C:153
Foam::Vector< scalar >
Foam::List< label >
Foam::valuePointPatchField::rmap
virtual void rmap(const pointPatchField< Type > &, const labelList &)
Reverse map the given PointPatchField onto.
Definition: valuePointPatchField.C:118
pointPatchFields.H
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::cos
dimensionedScalar cos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:265