wedgeFvPatchField.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-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 
26 #include "wedgeFvPatch.H"
27 #include "wedgeFvPatchField.H"
28 #include "transformField.H"
29 #include "symmTransform.H"
30 #include "diagTensor.H"
31 
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
35 template<class Type>
37 (
38  const fvPatch& p,
40 )
41 :
43 {}
44 
45 
46 template<class Type>
48 (
49  const wedgeFvPatchField<Type>& ptf,
50  const fvPatch& p,
52  const fvPatchFieldMapper& mapper
53 )
54 :
55  transformFvPatchField<Type>(ptf, p, iF, mapper)
56 {
57  if (!isType<wedgeFvPatch>(this->patch()))
58  {
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 (
72  const fvPatch& p,
74  const dictionary& dict
75 )
76 :
78 {
79  if (!isType<wedgeFvPatch>(p))
80  {
82  (
83  dict
84  ) << "\n patch type '" << p.type()
85  << "' not constraint type '" << typeName << "'"
86  << "\n for patch " << p.name()
87  << " of field " << this->internalField().name()
88  << " in file " << this->internalField().objectPath()
89  << exit(FatalIOError);
90  }
91 
92  evaluate();
93 }
94 
95 
96 template<class Type>
98 (
99  const wedgeFvPatchField<Type>& ptf
100 )
101 :
103 {}
104 
105 
106 template<class Type>
108 (
109  const wedgeFvPatchField<Type>& ptf,
111 )
112 :
114 {}
115 
116 
117 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118 
119 template<class Type>
121 {
122  const Field<Type> pif(this->patchInternalField());
123 
124  return
125  (
126  transform(refCast<const wedgeFvPatch>(this->patch()).cellT(), pif) - pif
127  )*(0.5*this->patch().deltaCoeffs());
128 }
129 
130 
131 template<class Type>
133 {
134  if (!this->updated())
135  {
136  this->updateCoeffs();
137  }
138 
140  (
141  transform
142  (
143  refCast<const wedgeFvPatch>(this->patch()).faceT(),
144  this->patchInternalField()
145  )
146  );
147 }
148 
149 
150 template<class Type>
153 {
154  const diagTensor diagT =
155  0.5*diag(I - refCast<const wedgeFvPatch>(this->patch()).cellT());
156 
157  const vector diagV(diagT.xx(), diagT.yy(), diagT.zz());
158 
159  return tmp<Field<Type>>
160  (
161  new Field<Type>
162  (
163  this->size(),
164  transformMask<Type>
165  (
166  pow
167  (
168  diagV,
170  ::type>::zero
171  )
172  )
173  )
174  );
175 }
176 
177 
178 // ************************************************************************* //
const Cmpt & xx() const
Definition: DiagTensorI.H:78
3D symmetric tensor transformation operations.
dictionary dict
const Cmpt & yy() const
Definition: DiagTensorI.H:84
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
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
virtual tmp< Field< Type > > snGradTransformDiag() const
Return face-gradient transform diagonal.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Traits class for primitives.
Definition: pTraits.H:50
const Cmpt & zz() const
Definition: DiagTensorI.H:90
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual tmp< Field< Type > > snGrad() const
Return gradient at boundary.
wedgeFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Spatial transformation functions for primitive fields.
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.
This boundary condition is similar to the cyclic condition, except that it is applied to 2-D geometri...
Foam::transformFvPatchField.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
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
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477
IOerror FatalIOError