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-2021 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 (
46  const fvPatch& p,
48  const dictionary& dict
49 )
50 :
52  symmetryPlanePatch_(refCast<const symmetryPlaneFvPatch>(p))
53 {
54  if (!isType<symmetryPlaneFvPatch>(p))
55  {
57  (
58  dict
59  ) << "\n patch type '" << p.type()
60  << "' not constraint type '" << typeName << "'"
61  << "\n for patch " << p.name()
62  << " of field " << this->internalField().name()
63  << " in file " << this->internalField().objectPath()
64  << exit(FatalIOError);
65  }
66 }
67 
68 
69 template<class Type>
71 (
73  const fvPatch& p,
75  const fvPatchFieldMapper& mapper
76 )
77 :
78  basicSymmetryFvPatchField<Type>(ptf, p, iF, mapper),
79  symmetryPlanePatch_(refCast<const symmetryPlaneFvPatch>(p))
80 {
81  if (!isType<symmetryPlaneFvPatch>(this->patch()))
82  {
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 (
98 )
99 :
101  symmetryPlanePatch_(ptf.symmetryPlanePatch_)
102 {}
103 
104 
105 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
106 
107 template<class Type>
110 {
111  vector nHat(symmetryPlanePatch_.n());
112 
113  const Field<Type> iF(this->patchInternalField());
114 
115  return
116  (transform(I - 2.0*sqr(nHat), iF) - iF)
117  *(this->patch().deltaCoeffs()/2.0);
118 }
119 
120 
121 template<class Type>
123 {
124  if (!this->updated())
125  {
126  this->updateCoeffs();
127  }
128 
129  vector nHat(symmetryPlanePatch_.n());
130 
131  const Field<Type> iF(this->patchInternalField());
132 
134  (
135  (iF + transform(I - 2.0*sqr(nHat), iF))/2.0
136  );
137 
139 }
140 
141 
142 template<class Type>
145 {
146  vector nHat(symmetryPlanePatch_.n());
147 
148  const vector diag
149  (
150  mag(nHat.component(vector::X)),
151  mag(nHat.component(vector::Y)),
152  mag(nHat.component(vector::Z))
153  );
154 
155  return tmp<Field<Type>>
156  (
157  new Field<Type>
158  (
159  this->size(),
160  transformMask<Type>
161  (
162  // pow<vector, pTraits<Type>::rank>(diag)
163  pow
164  (
165  diag,
167  ::type>::zero
168  )
169  )
170  )
171  );
172 }
173 
174 
175 // ************************************************************************* //
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
void evaluate(GeometricField< Type, PatchField, GeoMesh > &result, const Function1< Type > &func, const GeometricField< Type, PatchField, GeoMesh > &x)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
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:62
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.
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:335
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