symmetryPlaneFvPatchField.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) 2013-2018 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 
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Type>
33 (
34  const fvPatch& p,
36 )
37 :
39  symmetryPlanePatch_(refCast<const symmetryPlaneFvPatch>(p))
40 {}
41 
42 
43 template<class Type>
45 (
47  const fvPatch& p,
49  const fvPatchFieldMapper& mapper
50 )
51 :
52  basicSymmetryFvPatchField<Type>(ptf, p, iF, mapper),
53  symmetryPlanePatch_(refCast<const symmetryPlaneFvPatch>(p))
54 {
55  if (!isType<symmetryPlaneFvPatch>(this->patch()))
56  {
58  << "' not constraint type '" << typeName << "'"
59  << "\n for patch " << p.name()
60  << " of field " << this->internalField().name()
61  << " in file " << this->internalField().objectPath()
62  << exit(FatalIOError);
63  }
64 }
65 
66 
67 template<class Type>
69 (
70  const fvPatch& p,
72  const dictionary& dict
73 )
74 :
76  symmetryPlanePatch_(refCast<const symmetryPlaneFvPatch>(p))
77 {
78  if (!isType<symmetryPlaneFvPatch>(p))
79  {
81  (
82  dict
83  ) << "\n patch type '" << p.type()
84  << "' not constraint type '" << typeName << "'"
85  << "\n for patch " << p.name()
86  << " of field " << this->internalField().name()
87  << " in file " << this->internalField().objectPath()
88  << exit(FatalIOError);
89  }
90 }
91 
92 
93 template<class Type>
95 (
97 )
98 :
100  symmetryPlanePatch_(ptf.symmetryPlanePatch_)
101 {}
102 
103 
104 template<class Type>
106 (
109 )
110 :
112  symmetryPlanePatch_(ptf.symmetryPlanePatch_)
113 {}
114 
115 
116 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
117 
118 template<class Type>
121 {
122  vector nHat(symmetryPlanePatch_.n());
123 
124  const Field<Type> iF(this->patchInternalField());
125 
126  return
127  (transform(I - 2.0*sqr(nHat), iF) - iF)
128  *(this->patch().deltaCoeffs()/2.0);
129 }
130 
131 
132 template<class Type>
134 {
135  if (!this->updated())
136  {
137  this->updateCoeffs();
138  }
139 
140  vector nHat(symmetryPlanePatch_.n());
141 
142  const Field<Type> iF(this->patchInternalField());
143 
145  (
146  (iF + transform(I - 2.0*sqr(nHat), iF))/2.0
147  );
148 
150 }
151 
152 
153 template<class Type>
156 {
157  vector nHat(symmetryPlanePatch_.n());
158 
159  const vector diag
160  (
161  mag(nHat.component(vector::X)),
162  mag(nHat.component(vector::Y)),
163  mag(nHat.component(vector::Z))
164  );
165 
166  return tmp<Field<Type>>
167  (
168  new Field<Type>
169  (
170  this->size(),
171  transformMask<Type>
172  (
173  // pow<vector, pTraits<Type>::rank>(diag)
174  pow
175  (
176  diag,
178  ::type>::zero
179  )
180  )
181  )
182  );
183 }
184 
185 
186 // ************************************************************************* //
dictionary dict
virtual tmp< Field< Type > > snGradTransformDiag() const
Return face-gradient transform diagonal.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
commsTypes
Types of communications.
Definition: UPstream.H:64
dimensionedSymmTensor sqr(const dimensionedVector &dv)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Traits class for primitives.
Definition: pTraits.H:50
symmetryPlaneFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
static const Identity< scalar > I
Definition: Identity.H:93
Pre-declare SubField and related Field type.
Definition: Field.H:56
Foam::fvPatchFieldMapper.
virtual tmp< Field< Type > > snGrad() const
Return gradient at boundary.
Foam::transformFvPatchField.
This boundary condition enforces a symmetryPlane constraint.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
PtrList< volScalarField > & Y
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:49
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477
IOerror FatalIOError