wedgeFvPatchField.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "wedgeFvPatch.H"
29#include "wedgeFvPatchField.H"
30#include "transformField.H"
31#include "symmTransform.H"
32#include "diagTensor.H"
33
34
35// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36
37template<class Type>
39(
40 const fvPatch& p,
42)
43:
44 transformFvPatchField<Type>(p, iF)
45{}
46
47
48template<class Type>
50(
51 const wedgeFvPatchField<Type>& ptf,
52 const fvPatch& p,
54 const fvPatchFieldMapper& mapper
55)
56:
57 transformFvPatchField<Type>(ptf, p, iF, mapper)
58{
59 if (!isType<wedgeFvPatch>(this->patch()))
60 {
62 << "' not constraint type '" << typeName << "'"
63 << "\n for patch " << p.name()
64 << " of field " << this->internalField().name()
65 << " in file " << this->internalField().objectPath()
66 << exit(FatalError);
67 }
68}
69
70
71template<class Type>
73(
74 const fvPatch& p,
76 const dictionary& dict
77)
78:
79 transformFvPatchField<Type>(p, iF, dict)
80{
81 if (!isType<wedgeFvPatch>(p))
82 {
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()
90 }
91
92 evaluate();
93}
94
95
96template<class Type>
98(
100)
101:
102 transformFvPatchField<Type>(ptf)
103{}
104
105
106template<class Type>
108(
109 const wedgeFvPatchField<Type>& ptf,
111)
112:
113 transformFvPatchField<Type>(ptf, iF)
114{}
115
116
117// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118
119template<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
131template<class Type>
133{
134 if (!this->updated())
135 {
136 this->updateCoeffs();
137 }
138
140 (
142 (
143 refCast<const wedgeFvPatch>(this->patch()).faceT(),
144 this->patchInternalField()
145 )
146 );
147}
148
149
150template<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:76
const Cmpt & zz() const
Definition: DiagTensorI.H:88
const Cmpt & yy() const
Definition: DiagTensorI.H:82
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static const char *const typeName
Typename for Field.
Definition: FieldBase.H:59
Generic templated field type.
Definition: Field.H:82
void evaluate()
Evaluate boundary conditions.
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
commsTypes
Types of communications.
Definition: UPstream.H:67
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A FieldMapper for finite-volume patch fields.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:82
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:368
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:362
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
virtual const word & name() const
Return name.
Definition: fvPatch.H:173
A traits class, which is primarily used for primitives.
Definition: pTraits.H:59
A class for managing temporary objects.
Definition: tmp.H:65
Foam::transformFvPatchField.
This boundary condition is similar to the cyclic condition, except that it is applied to 2-D geometri...
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
virtual tmp< Field< Type > > snGrad() const
Return gradient at boundary.
virtual tmp< Field< Type > > snGradTransformDiag() const
Return face-gradient transform diagonal.
volScalarField & p
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Definition: dimensionSet.C:536
static const Identity< scalar > I
Definition: Identity.H:94
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
IOerror FatalIOError
error FatalError
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dictionary dict
3D symmetric tensor transformation operations.
Spatial transformation functions for primitive fields.