symmetryPlanePointPatchField.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-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 
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Type>
32 (
33  const pointPatch& p,
35 )
36 :
38  symmetryPlanePatch_(refCast<const symmetryPlanePointPatch>(p))
39 {}
40 
41 
42 template<class Type>
44 (
45  const pointPatch& p,
47  const dictionary& dict
48 )
49 :
51  symmetryPlanePatch_(refCast<const symmetryPlanePointPatch>(p))
52 {
53  if (!isType<symmetryPlanePointPatch>(p))
54  {
56  (
57  dict
58  ) << "patch " << this->patch().index() << " not symmetry type. "
59  << "Patch type = " << p.type()
60  << exit(FatalIOError);
61  }
62 }
63 
64 
65 template<class Type>
67 (
69  const pointPatch& p,
71  const fieldMapper& mapper
72 )
73 :
74  basicSymmetryPointPatchField<Type>(ptf, p, iF, mapper),
75  symmetryPlanePatch_(refCast<const symmetryPlanePointPatch>(p))
76 {
77  if (!isType<symmetryPlanePointPatch>(this->patch()))
78  {
80  << "Field type does not correspond to patch type for patch "
81  << this->patch().index() << "." << endl
82  << "Field type: " << typeName << endl
83  << "Patch type: " << this->patch().type()
84  << exit(FatalError);
85  }
86 }
87 
88 
89 template<class Type>
91 (
94 )
95 :
96  basicSymmetryPointPatchField<Type>(ptf, iF),
97  symmetryPlanePatch_(ptf.symmetryPlanePatch_)
98 {}
99 
100 
101 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
102 
103 template<class Type>
105 (
106  const Pstream::commsTypes
107 )
108 {
109  vector nHat = symmetryPlanePatch_.n();
110 
111  tmp<Field<Type>> tvalues =
112  (
113  (
114  this->patchInternalField()
115  + transform(I - 2.0*sqr(nHat), this->patchInternalField())
116  )/2.0
117  );
118 
119  // Get internal field to insert values into
120  Field<Type>& iF = const_cast<Field<Type>&>(this->primitiveField());
121 
122  this->setInternalField(iF, tvalues());
123 }
124 
125 
126 // ************************************************************************* //
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
commsTypes
Types of communications.
Definition: UPstream.H:65
A Symmetry boundary condition for pointField.
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
const pointPatch & patch() const
Return patch.
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:57
virtual label index() const =0
Return the index of this patch in the pointBoundaryMesh.
A symmetry-plane boundary condition for pointField.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Update the patch field.
symmetryPlanePointPatchField(const pointPatch &, const DimensionedField< Type, pointMesh > &)
Construct from patch and internal field.
A class for managing temporary objects.
Definition: tmp.H:55
#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
dimensionedSymmTensor sqr(const dimensionedVector &dv)
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:129
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
static const Identity< scalar > I
Definition: Identity.H:93
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
IOerror FatalIOError
error FatalError
dictionary dict
volScalarField & p