directionMixedFvPatchField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2024 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Type>
33 (
34  const fvPatch& p,
36 )
37 :
38  transformFvPatchField<Type>(p, iF),
39  refValue_(p.size()),
40  refGrad_(p.size()),
41  valueFraction_(p.size())
42 {}
43 
44 
45 template<class Type>
47 (
48  const fvPatch& p,
50  const dictionary& dict
51 )
52 :
53  transformFvPatchField<Type>(p, iF, dict),
54  refValue_("refValue", iF.dimensions(), dict, p.size()),
55  refGrad_("refGradient", iF.dimensions()/dimLength, dict, p.size()),
56  valueFraction_("valueFraction", unitFraction, dict, p.size())
57 {
58  evaluate();
59 }
60 
61 
62 template<class Type>
64 (
66  const fvPatch& p,
68  const fieldMapper& mapper
69 )
70 :
71  transformFvPatchField<Type>(ptf, p, iF, mapper),
72  refValue_(mapper(ptf.refValue_)),
73  refGrad_(mapper(ptf.refGrad_)),
74  valueFraction_(mapper(ptf.valueFraction_))
75 {}
76 
77 
78 template<class Type>
80 (
83 )
84 :
85  transformFvPatchField<Type>(ptf, iF),
86  refValue_(ptf.refValue_),
87  refGrad_(ptf.refGrad_),
88  valueFraction_(ptf.valueFraction_)
89 {}
90 
91 
92 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
93 
94 template<class Type>
96 (
97  const fvPatchField<Type>& ptf,
98  const fieldMapper& mapper
99 )
100 {
102 
103  const directionMixedFvPatchField<Type>& dmptf =
104  refCast<const directionMixedFvPatchField<Type>>(ptf);
105 
106  mapper(refValue_, dmptf.refValue_);
107  mapper(refGrad_, dmptf.refGrad_);
108  mapper(valueFraction_, dmptf.valueFraction_);
109 }
110 
111 
112 template<class Type>
114 (
115  const fvPatchField<Type>& ptf
116 )
117 {
119 
120  const directionMixedFvPatchField<Type>& dmptf =
121  refCast<const directionMixedFvPatchField<Type>>(ptf);
122 
123  refValue_.reset(dmptf.refValue_);
124  refGrad_.reset(dmptf.refGrad_);
125  valueFraction_.reset(dmptf.valueFraction_);
126 }
127 
128 
129 template<class Type>
132 {
133  const Field<Type> pif(this->patchInternalField());
134 
135  tmp<Field<Type>> normalValue = transform(valueFraction_, refValue_);
136 
137  tmp<Field<Type>> gradValue = pif + refGrad_/this->patch().deltaCoeffs();
138 
139  tmp<Field<Type>> transformGradValue =
140  transform(I - valueFraction_, gradValue);
141 
142  return
143  (normalValue + transformGradValue - pif)*
144  this->patch().deltaCoeffs();
145 }
146 
147 
148 template<class Type>
150 {
151  if (!this->updated())
152  {
153  this->updateCoeffs();
154  }
155 
156  tmp<Field<Type>> normalValue = transform(valueFraction_, refValue_);
157 
158  tmp<Field<Type>> gradValue =
159  this->patchInternalField() + refGrad_/this->patch().deltaCoeffs();
160 
161  tmp<Field<Type>> transformGradValue =
162  transform(I - valueFraction_, gradValue);
163 
164  Field<Type>::operator=(normalValue + transformGradValue);
165 
167 }
168 
169 
170 template<class Type>
173 {
174  vectorField diag(valueFraction_.size());
175 
176  diag.replace
177  (
178  vector::X,
179  sqrt(mag(valueFraction_.component(symmTensor::XX)))
180  );
181  diag.replace
182  (
183  vector::Y,
184  sqrt(mag(valueFraction_.component(symmTensor::YY)))
185  );
186  diag.replace
187  (
188  vector::Z,
189  sqrt(mag(valueFraction_.component(symmTensor::ZZ)))
190  );
191 
192  return transformFieldMask<Type>(pow<vector, pTraits<Type>::rank>(diag));
193 }
194 
195 
196 template<class Type>
198 {
200  writeEntry(os, "refValue", refValue_);
201  writeEntry(os, "refGradient", refGrad_);
202  writeEntry(os, "valueFraction", valueFraction_);
203  writeEntry(os, "value", *this);
204 }
205 
206 
207 // ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Pre-declare SubField and related Field type.
Definition: Field.H:83
void operator=(const Field< Type > &)
Definition: Field.C:550
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
commsTypes
Types of communications.
Definition: UPstream.H:65
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Base class for direction-mixed boundary conditions.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
virtual void write(Ostream &) const
Write.
virtual tmp< Field< Type > > snGrad() const
Return gradient at boundary.
virtual void reset(const fvPatchField< Type > &)
Reset the fvPatchField to the given fvPatchField.
virtual tmp< Field< Type > > snGradTransformDiag() const
Return face-gradient transform diagonal.
virtual void map(const fvPatchField< Type > &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
directionMixedFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Abstract base class for field mapping.
Definition: fieldMapper.H:48
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field, sets Updated to false.
Definition: fvPatchField.C:209
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:229
virtual void reset(const fvPatchField< Type > &)
Reset the fvPatchField to the given fvPatchField.
Definition: fvPatchField.C:195
virtual void map(const fvPatchField< Type > &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:185
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
Traits class for primitives.
Definition: pTraits.H:53
A class for managing temporary objects.
Definition: tmp.H:55
Foam::transformFvPatchField.
const HashTable< dimensionSet > & dimensions()
Get the table of dimension sets.
Definition: dimensionSets.C:96
const dimensionSet dimLength
static const Identity< scalar > I
Definition: Identity.H:93
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionedScalar sqrt(const dimensionedScalar &ds)
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
dimensioned< scalar > mag(const dimensioned< Type > &)
const unitConversion unitFraction
dictionary dict
volScalarField & p
Spatial transformation functions for symmTensor fields.