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-2023 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 :
42  transformFvPatchField<Type>(p, iF)
43 {}
44 
45 
46 template<class Type>
48 (
49  const fvPatch& p,
51  const dictionary& dict
52 )
53 :
54  transformFvPatchField<Type>(p, iF, dict)
55 {
56  if (!isType<wedgeFvPatch>(p))
57  {
59  (
60  dict
61  ) << "\n patch type '" << p.type()
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(FatalIOError);
67  }
68 
69  evaluate();
70 }
71 
72 
73 template<class Type>
75 (
76  const wedgeFvPatchField<Type>& ptf,
77  const fvPatch& p,
79  const fieldMapper& mapper
80 )
81 :
82  transformFvPatchField<Type>(ptf, p, iF, mapper)
83 {
84  if (!isType<wedgeFvPatch>(this->patch()))
85  {
87  << "' not constraint type '" << typeName << "'"
88  << "\n for patch " << p.name()
89  << " of field " << this->internalField().name()
90  << " in file " << this->internalField().objectPath()
91  << exit(FatalIOError);
92  }
93 }
94 
95 
96 template<class Type>
98 (
99  const wedgeFvPatchField<Type>& ptf,
101 )
102 :
103  transformFvPatchField<Type>(ptf, iF)
104 {}
105 
106 
107 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
108 
109 template<class Type>
111 {
112  const Field<Type> pif(this->patchInternalField());
113 
114  return
115  (
116  transform(refCast<const wedgeFvPatch>(this->patch()).cellT(), pif) - pif
117  )*(0.5*this->patch().deltaCoeffs());
118 }
119 
120 
121 template<class Type>
123 {
124  if (!this->updated())
125  {
126  this->updateCoeffs();
127  }
128 
130  (
131  transform
132  (
133  refCast<const wedgeFvPatch>(this->patch()).faceT(),
134  this->patchInternalField()
135  )
136  );
137 }
138 
139 
140 template<class Type>
143 {
144  const vector diagV
145  (
146  0.5*diag(I - refCast<const wedgeFvPatch>(this->patch()).cellT())
147  );
148 
149  return tmp<Field<Type>>
150  (
151  new Field<Type>
152  (
153  this->size(),
154  transformMask<Type>
155  (
156  pow
157  (
158  diagV,
160  ::type>::zero
161  )
162  )
163  )
164  );
165 }
166 
167 
168 // ************************************************************************* //
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
static const char *const typeName
Definition: Field.H:106
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
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
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:362
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:356
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
virtual const word & name() const
Return name.
Definition: fvPatch.H:126
Traits class for primitives.
Definition: pTraits.H:53
A class for managing temporary objects.
Definition: tmp.H:55
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.
wedgeFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
virtual tmp< Field< Type > > snGrad() const
Return gradient at boundary.
virtual tmp< Field< Type > > snGradTransformDiag() const
Return face-gradient transform diagonal.
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:50
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:346
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const Identity< scalar > I
Definition: Identity.H:93
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
IOerror FatalIOError
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict
volScalarField & p
3D symmetric tensor transformation operations.
Spatial transformation functions for primitive fields.