mixedFixedValueSlipFvPatchField.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-2012 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 "symmTransformField.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
36 template<class Type>
38 (
39  const fvPatch& p,
40  const DimensionedField<Type, volMesh>& iF
41 )
42 :
43  transformFvPatchField<Type>(p, iF),
44  refValue_(p.size()),
45  valueFraction_(p.size(), 1.0)
46 {}
47 
48 
49 template<class Type>
50 mixedFixedValueSlipFvPatchField<Type>::mixedFixedValueSlipFvPatchField
51 (
52  const mixedFixedValueSlipFvPatchField<Type>& ptf,
53  const fvPatch& p,
54  const DimensionedField<Type, volMesh>& iF,
55  const fvPatchFieldMapper& mapper
56 )
57 :
58  transformFvPatchField<Type>(ptf, p, iF, mapper),
59  refValue_(ptf.refValue_, mapper),
60  valueFraction_(ptf.valueFraction_, mapper)
61 {}
62 
63 
64 template<class Type>
65 mixedFixedValueSlipFvPatchField<Type>::mixedFixedValueSlipFvPatchField
66 (
67  const fvPatch& p,
68  const DimensionedField<Type, volMesh>& iF,
69  const dictionary& dict
70 )
71 :
72  transformFvPatchField<Type>(p, iF),
73  refValue_("refValue", dict, p.size()),
74  valueFraction_("valueFraction", dict, p.size())
75 {}
76 
77 
78 template<class Type>
79 mixedFixedValueSlipFvPatchField<Type>::mixedFixedValueSlipFvPatchField
80 (
81  const mixedFixedValueSlipFvPatchField<Type>& ptf
82 )
83 :
84  transformFvPatchField<Type>(ptf),
85  refValue_(ptf.refValue_),
86  valueFraction_(ptf.valueFraction_)
87 {}
88 
89 template<class Type>
90 mixedFixedValueSlipFvPatchField<Type>::mixedFixedValueSlipFvPatchField
91 (
92  const mixedFixedValueSlipFvPatchField<Type>& ptf,
93  const DimensionedField<Type, volMesh>& iF
94 )
95 :
96  transformFvPatchField<Type>(ptf, iF),
97  refValue_(ptf.refValue_),
98  valueFraction_(ptf.valueFraction_)
99 {}
100 
101 
102 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
103 
104 // Map from self
105 template<class Type>
106 void mixedFixedValueSlipFvPatchField<Type>::autoMap
107 (
108  const fvPatchFieldMapper& m
109 )
110 {
111  Field<Type>::autoMap(m);
112  refValue_.autoMap(m);
113  valueFraction_.autoMap(m);
114 }
115 
116 
117 // Reverse-map the given fvPatchField onto this fvPatchField
118 template<class Type>
119 void mixedFixedValueSlipFvPatchField<Type>::rmap
120 (
121  const fvPatchField<Type>& ptf,
122  const labelList& addr
123 )
124 {
125  transformFvPatchField<Type>::rmap(ptf, addr);
126 
127  const mixedFixedValueSlipFvPatchField<Type>& dmptf =
128  refCast<const mixedFixedValueSlipFvPatchField<Type> >(ptf);
129 
130  refValue_.rmap(dmptf.refValue_, addr);
131  valueFraction_.rmap(dmptf.valueFraction_, addr);
132 }
133 
134 
135 // Return gradient at boundary
136 template<class Type>
137 tmp<Field<Type> > mixedFixedValueSlipFvPatchField<Type>::snGrad() const
138 {
139  tmp<vectorField> nHat = this->patch().nf();
140  Field<Type> pif(this->patchInternalField());
141 
142  return
143  (
144  valueFraction_*refValue_
145  + (1.0 - valueFraction_)*transform(I - sqr(nHat), pif) - pif
146  )*this->patch().deltaCoeffs();
147 }
148 
149 
150 // Evaluate the field on the patch
151 template<class Type>
152 void mixedFixedValueSlipFvPatchField<Type>::evaluate(const Pstream::commsTypes)
153 {
154  if (!this->updated())
155  {
156  this->updateCoeffs();
157  }
158 
159  vectorField nHat(this->patch().nf());
160 
161  Field<Type>::operator=
162  (
163  valueFraction_*refValue_
164  +
165  (1.0 - valueFraction_)
166  *transform(I - nHat*nHat, this->patchInternalField())
167  );
168 
169  transformFvPatchField<Type>::evaluate();
170 }
171 
172 
173 // Return defining fields
174 template<class Type>
175 tmp<Field<Type> >
176 mixedFixedValueSlipFvPatchField<Type>::snGradTransformDiag() const
177 {
178  vectorField nHat(this->patch().nf());
179  vectorField diag(nHat.size());
180 
181  diag.replace(vector::X, mag(nHat.component(vector::X)));
182  diag.replace(vector::Y, mag(nHat.component(vector::Y)));
183  diag.replace(vector::Z, mag(nHat.component(vector::Z)));
184 
185  return
186  valueFraction_*Type(pTraits<Type>::one)
187  + (1.0 - valueFraction_)
188  *transformFieldMask<Type>(pow<vector, pTraits<Type>::rank>(diag));
189 }
190 
191 
192 // Write
193 template<class Type>
194 void mixedFixedValueSlipFvPatchField<Type>::write(Ostream& os) const
195 {
197  refValue_.writeEntry("refValue", os);
198  valueFraction_.writeEntry("valueFraction", os);
199 }
200 
201 
202 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
203 
204 } // End namespace Foam
205 
206 // ************************************************************************* //
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > mag(const dimensioned< Type > &)
static const sphericalTensor I(1)
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
Namespace for OpenFOAM.
runTime write()
dictionary dict
volVectorField vectorField(fieldObject, mesh)
volScalarField & p
Definition: createFields.H:51
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:465
PtrList< volScalarField > & Y
Definition: createFields.H:36
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
List< label > labelList
A List of labels.
Definition: labelList.H:56
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
mixedFixedValueSlipFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)