mixedFvPatchField.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-2016 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 
26 #include "mixedFvPatchField.H"
27 
28 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29 
30 template<class Type>
32 (
33  const fvPatch& p,
35 )
36 :
37  fvPatchField<Type>(p, iF),
38  refValue_(p.size()),
39  refGrad_(p.size()),
40  valueFraction_(p.size())
41 {}
42 
43 
44 template<class Type>
46 (
47  const fvPatch& p,
49  const dictionary& dict
50 )
51 :
53  refValue_("refValue", dict, p.size()),
54  refGrad_("refGradient", dict, p.size()),
55  valueFraction_("valueFraction", dict, p.size())
56 {
57  evaluate();
58 }
59 
60 
61 template<class Type>
63 (
64  const mixedFvPatchField<Type>& ptf,
65  const fvPatch& p,
67  const fvPatchFieldMapper& mapper
68 )
69 :
70  fvPatchField<Type>(ptf, p, iF, mapper),
71  refValue_(ptf.refValue_, mapper),
72  refGrad_(ptf.refGrad_, mapper),
73  valueFraction_(ptf.valueFraction_, mapper)
74 {
75  if (notNull(iF) && mapper.hasUnmapped())
76  {
78  << "On field " << iF.name() << " patch " << p.name()
79  << " patchField " << this->type()
80  << " : mapper does not map all values." << nl
81  << " To avoid this warning fully specify the mapping in derived"
82  << " patch fields." << endl;
83  }
84 }
85 
86 
87 template<class Type>
89 (
90  const mixedFvPatchField<Type>& ptf
91 )
92 :
93  fvPatchField<Type>(ptf),
94  refValue_(ptf.refValue_),
95  refGrad_(ptf.refGrad_),
96  valueFraction_(ptf.valueFraction_)
97 {}
98 
99 
100 template<class Type>
102 (
103  const mixedFvPatchField<Type>& ptf,
105 )
106 :
107  fvPatchField<Type>(ptf, iF),
108  refValue_(ptf.refValue_),
109  refGrad_(ptf.refGrad_),
110  valueFraction_(ptf.valueFraction_)
111 {}
112 
113 
114 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
115 
116 template<class Type>
118 (
119  const fvPatchFieldMapper& m
120 )
121 {
123  refValue_.autoMap(m);
124  refGrad_.autoMap(m);
125  valueFraction_.autoMap(m);
126 }
127 
128 
129 template<class Type>
131 (
132  const fvPatchField<Type>& ptf,
133  const labelList& addr
134 )
135 {
136  fvPatchField<Type>::rmap(ptf, addr);
137 
138  const mixedFvPatchField<Type>& mptf =
139  refCast<const mixedFvPatchField<Type>>(ptf);
140 
141  refValue_.rmap(mptf.refValue_, addr);
142  refGrad_.rmap(mptf.refGrad_, addr);
143  valueFraction_.rmap(mptf.valueFraction_, addr);
144 }
145 
146 
147 template<class Type>
149 {
150  if (!this->updated())
151  {
152  this->updateCoeffs();
153  }
154 
156  (
157  valueFraction_*refValue_
158  +
159  (1.0 - valueFraction_)*
160  (
161  this->patchInternalField()
162  + refGrad_/this->patch().deltaCoeffs()
163  )
164  );
165 
167 }
168 
169 
170 template<class Type>
173 {
174  return
175  valueFraction_
176  *(refValue_ - this->patchInternalField())
177  *this->patch().deltaCoeffs()
178  + (1.0 - valueFraction_)*refGrad_;
179 }
180 
181 
182 template<class Type>
185 (
186  const tmp<scalarField>&
187 ) const
188 {
189  return Type(pTraits<Type>::one)*(1.0 - valueFraction_);
190 
191 }
192 
193 
194 template<class Type>
197 (
198  const tmp<scalarField>&
199 ) const
200 {
201  return
202  valueFraction_*refValue_
203  + (1.0 - valueFraction_)*refGrad_/this->patch().deltaCoeffs();
204 }
205 
206 
207 template<class Type>
210 {
211  return -Type(pTraits<Type>::one)*valueFraction_*this->patch().deltaCoeffs();
212 }
213 
214 
215 template<class Type>
218 {
219  return
220  valueFraction_*this->patch().deltaCoeffs()*refValue_
221  + (1.0 - valueFraction_)*refGrad_;
222 }
223 
224 
225 template<class Type>
227 {
229  refValue_.writeEntry("refValue", os);
230  refGrad_.writeEntry("refGradient", os);
231  valueFraction_.writeEntry("valueFraction", os);
232  this->writeEntry("value", os);
233 }
234 
235 
236 // ************************************************************************* //
dictionary dict
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
commsTypes
Types of communications.
Definition: UPstream.H:64
virtual tmp< Field< Type > > snGrad() const
Return gradient at boundary.
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual tmp< Field< Type > > gradientInternalCoeffs() const
Return the matrix diagonal coefficients corresponding to the.
Traits class for primitives.
Definition: pTraits.H:50
This boundary condition provides a base class for &#39;mixed&#39; type boundary conditions, i.e. conditions that mix fixed value and patch-normal gradient conditions.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::blocking)
Evaluate the patch field.
Pre-declare SubField and related Field type.
Definition: Field.H:57
Foam::fvPatchFieldMapper.
virtual tmp< Field< Type > > valueInternalCoeffs(const tmp< scalarField > &) const
Return the matrix diagonal coefficients corresponding to the.
virtual bool hasUnmapped() const =0
Are there unmapped values? I.e. do all size() elements get.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:262
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual tmp< Field< Type > > gradientBoundaryCoeffs() const
Return the matrix source coefficients corresponding to the.
bool notNull(const T &t)
Return true if t is not a reference to the nullObject of type T.
Definition: nullObjectI.H:46
virtual void write(Ostream &) const
Write.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
#define WarningInFunction
Report a warning using Foam::Warning.
mixedFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual tmp< Field< Type > > valueBoundaryCoeffs(const tmp< scalarField > &) const
Return the matrix source coefficients corresponding to the.
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:54
runTime write()
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.