oscillatingFixedValueFvPatchField.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 | Copyright (C) 2011-2013 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 "mathematicalConstants.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 
34 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
35 
36 template<class Type>
37 scalar oscillatingFixedValueFvPatchField<Type>::currentScale() const
38 {
39  const scalar t = this->db().time().timeOutputValue();
40  const scalar a = amplitude_->value(t);
41  const scalar f = frequency_->value(t);
42 
43  return 1.0 + a*sin(constant::mathematical::twoPi*f*t);
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
49 template<class Type>
51 (
52  const fvPatch& p,
54 )
55 :
57  refValue_(p.size()),
58  offset_(pTraits<Type>::zero),
59  amplitude_(),
60  frequency_(),
61  curTimeIndex_(-1)
62 {}
63 
64 
65 template<class Type>
67 (
69  const fvPatch& p,
71  const fvPatchFieldMapper& mapper
72 )
73 :
74  fixedValueFvPatchField<Type>(ptf, p, iF, mapper),
75  refValue_(ptf.refValue_, mapper),
76  offset_(ptf.offset_),
77  amplitude_(ptf.amplitude_().clone().ptr()),
78  frequency_(ptf.frequency_().clone().ptr()),
79  curTimeIndex_(-1)
80 {}
81 
82 
83 template<class Type>
85 (
86  const fvPatch& p,
88  const dictionary& dict
89 )
90 :
92  refValue_("refValue", dict, p.size()),
93  offset_(dict.lookupOrDefault<Type>("offset", pTraits<Type>::zero)),
94  amplitude_(DataEntry<scalar>::New("amplitude", dict)),
95  frequency_(DataEntry<scalar>::New("frequency", dict)),
96  curTimeIndex_(-1)
97 {
98  if (dict.found("value"))
99  {
101  (
102  Field<Type>("value", dict, p.size())
103  );
104  }
105  else
106  {
108  (
109  refValue_*currentScale()
110  + offset_
111  );
112  }
113 }
114 
115 
116 template<class Type>
118 (
120 )
121 :
123  refValue_(ptf.refValue_),
124  offset_(ptf.offset_),
125  amplitude_(ptf.amplitude_().clone().ptr()),
126  frequency_(ptf.frequency_().clone().ptr()),
127  curTimeIndex_(-1)
128 {}
129 
130 
131 template<class Type>
133 (
136 )
137 :
139  refValue_(ptf.refValue_),
140  offset_(ptf.offset_),
141  amplitude_(ptf.amplitude_().clone().ptr()),
142  frequency_(ptf.frequency_().clone().ptr()),
143  curTimeIndex_(-1)
144 {}
145 
146 
147 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
148 
149 template<class Type>
151 (
152  const fvPatchFieldMapper& m
153 )
154 {
156  refValue_.autoMap(m);
157 }
158 
159 
160 template<class Type>
162 (
163  const fvPatchField<Type>& ptf,
164  const labelList& addr
165 )
166 {
168 
170  refCast<const oscillatingFixedValueFvPatchField<Type> >(ptf);
171 
172  refValue_.rmap(tiptf.refValue_, addr);
173 }
174 
175 
176 template<class Type>
178 {
179  if (this->updated())
180  {
181  return;
182  }
183 
184  if (curTimeIndex_ != this->db().time().timeIndex())
185  {
187  (
188  refValue_*currentScale()
189  + offset_
190  );
191 
192  curTimeIndex_ = this->db().time().timeIndex();
193  }
194 
196 }
197 
198 
199 template<class Type>
201 {
203  refValue_.writeEntry("refValue", os);
204  os.writeKeyword("offset") << offset_ << token::END_STATEMENT << nl;
205  amplitude_->writeData(os);
206  frequency_->writeData(os);
207 }
208 
209 
210 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
211 
212 } // End namespace Foam
213 
214 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:301
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
oscillatingFixedValueFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
label timeIndex
Definition: getTimeIndex.H:4
void size(const label)
Override size to be inconsistent with allocated storage.
Foam::fvPatchFieldMapper.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Namespace for OpenFOAM.
dictionary dict
static const char nl
Definition: Ostream.H:260
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:291
volScalarField & p
Definition: createFields.H:51
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: fvPatchField.C:230
This boundary condition provides an oscillating condition in terms of amplitude and frequency...
Traits class for primitives.
Definition: pTraits.H:50
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
virtual void write(Ostream &) const
Write.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: DataEntry.H:52
const scalar twoPi(2 *pi)
virtual tmp< fvPatchField< Type > > clone() const
Construct and return a clone.
dimensionedScalar sin(const dimensionedScalar &ds)