FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 
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 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
102 
103 // Map from self
104 template<class Type>
105 void mixedFixedValueSlipFvPatchField<Type>::autoMap
106 (
107  const fvPatchFieldMapper& m
108 )
109 {
110  Field<Type>::autoMap(m);
111  refValue_.autoMap(m);
112  valueFraction_.autoMap(m);
113 }
114 
115 
116 // Reverse-map the given fvPatchField onto this fvPatchField
117 template<class Type>
118 void mixedFixedValueSlipFvPatchField<Type>::rmap
119 (
120  const fvPatchField<Type>& ptf,
121  const labelList& addr
122 )
123 {
124  transformFvPatchField<Type>::rmap(ptf, addr);
125 
126  const mixedFixedValueSlipFvPatchField<Type>& dmptf =
127  refCast<const mixedFixedValueSlipFvPatchField<Type> >(ptf);
128 
129  refValue_.rmap(dmptf.refValue_, addr);
130  valueFraction_.rmap(dmptf.valueFraction_, addr);
131 }
132 
133 
134 // Return gradient at boundary
135 template<class Type>
136 tmp<Field<Type> > mixedFixedValueSlipFvPatchField<Type>::snGrad() const
137 {
138  vectorField nHat = this->patch().nf();
139  Field<Type> pif = this->patchInternalField();
140 
141  return
142  (
143  valueFraction_*refValue_
144  + (1.0 - valueFraction_)*transform(I - sqr(nHat), pif) - pif
145  )*this->patch().deltaCoeffs();
146 }
147 
148 
149 // Evaluate the field on the patch
150 template<class Type>
151 void mixedFixedValueSlipFvPatchField<Type>::evaluate(const Pstream::commsTypes)
152 {
153  if (!this->updated())
154  {
155  this->updateCoeffs();
156  }
157 
158  vectorField nHat = this->patch().nf();
159 
160  Field<Type>::operator=
161  (
162  valueFraction_*refValue_
163  +
164  (1.0 - valueFraction_)
165  *transform(I - nHat*nHat, this->patchInternalField())
166  );
167 
168  transformFvPatchField<Type>::evaluate();
169 }
170 
171 
172 // Return defining fields
173 template<class Type>
174 tmp<Field<Type> > mixedFixedValueSlipFvPatchField<Type>::snGradTransformDiag() const
175 {
176  vectorField nHat = this->patch().nf();
177  vectorField diag(nHat.size());
178 
179  diag.replace(vector::X, mag(nHat.component(vector::X)));
180  diag.replace(vector::Y, mag(nHat.component(vector::Y)));
181  diag.replace(vector::Z, mag(nHat.component(vector::Z)));
182 
183  return
184  valueFraction_*Type(pTraits<Type>::one)
185  + (1.0 - valueFraction_)*transformFieldMask<Type>(pow<vector, pTraits<Type>::rank>(diag));
186 }
187 
188 
189 // Write
190 template<class Type>
191 void mixedFixedValueSlipFvPatchField<Type>::write(Ostream& os) const
192 {
193  transformFvPatchField<Type>::write(os);
194  refValue_.writeEntry("refValue", os);
195  valueFraction_.writeEntry("valueFraction", os);
196 }
197 
198 
199 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
200 
201 } // End namespace Foam
202 
203 // ************************ vim: set sw=4 sts=4 et: ************************ //