patchMeanVelocityForce.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) 2015-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 
26 #include "patchMeanVelocityForce.H"
27 #include "volFields.H"
30 
31 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace fv
36 {
37  defineTypeNameAndDebug(patchMeanVelocityForce, 0);
38 
40  (
41  fvConstraint,
42  patchMeanVelocityForce,
43  dictionary
44  );
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 void Foam::fv::patchMeanVelocityForce::readCoeffs()
52 {
53  patch_ = coeffs().lookup<word>("patch");
54 
55  if (mesh().boundaryMesh().findPatchID(patch_) < 0)
56  {
58  << "Cannot find patch " << patch_
59  << exit(FatalError);
60  }
61 }
62 
63 
64 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
65 
67 (
68  const word& sourceName,
69  const word& modelType,
70  const dictionary& dict,
71  const fvMesh& mesh
72 )
73 :
74  meanVelocityForce(sourceName, modelType, dict, mesh),
75  patch_(word::null)
76 {
77  readCoeffs();
78 }
79 
80 
81 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
82 
83 Foam::scalar Foam::fv::patchMeanVelocityForce::magUbarAve
84 (
85  const volVectorField& U
86 ) const
87 {
88  const label patchi = mesh().boundaryMesh().findPatchID(patch_);
89 
90  scalar sumA = sum(mesh().boundary()[patchi].magSf());
91  scalar sumAmagU =
92  sum
93  (
94  mesh().boundary()[patchi].magSf()
95  *(normalised(Ubar()) & U.boundaryField()[patchi])
96  );
97 
98  // If the mean velocity force is applied to a cyclic patch
99  // for parallel runs include contributions from processorCyclic patches
100  // generated from the decomposition of the cyclic patch
101  const polyBoundaryMesh& patches = mesh().boundaryMesh();
102  if (Pstream::parRun() && isA<cyclicPolyPatch>(patches[patchi]))
103  {
104  const labelList processorCyclicPatches
105  (
106  processorCyclicPolyPatch::patchIDs(patch_, patches)
107  );
108 
109  forAll(processorCyclicPatches, pcpi)
110  {
111  const label patchi = processorCyclicPatches[pcpi];
112 
113  sumA += sum(mesh().boundary()[patchi].magSf());
114  sumAmagU +=
115  sum
116  (
117  mesh().boundary()[patchi].magSf()
118  *(normalised(Ubar()) & U.boundaryField()[patchi])
119  );
120  }
121  }
122 
123  mesh().reduce(sumA, sumOp<scalar>());
124  mesh().reduce(sumAmagU, sumOp<scalar>());
125 
126  return sumAmagU/sumA;
127 }
128 
129 
131 {
132  if (meanVelocityForce::read(dict))
133  {
134  readCoeffs();
135  return true;
136  }
137  else
138  {
139  return false;
140  }
141 }
142 
143 
144 // ************************************************************************* //
patchMeanVelocityForce(const word &sourceName, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from explicit source name and mesh.
defineTypeNameAndDebug(fixedTemperatureConstraint, 0)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
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
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual bool read(const dictionary &dict)
Read source dictionary.
patches[0]
Macros for easy insertion into run-time selection tables.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
dynamicFvMesh & mesh
Form normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
Definition: VectorSpaceI.H:413
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
static const word null
An empty word.
Definition: word.H:77
faceListList boundary(nPatches)
Foam::polyBoundaryMesh.
addToRunTimeSelectionTable(fvConstraint, fixedTemperatureConstraint, dictionary)
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
label patchi
dimensionedVector Ubar("Ubar", dimVelocity, laminarTransport)
virtual bool read(const dictionary &dict)
Read dictionary.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
static labelList patchIDs(const word &cyclicPolyPatchName, const polyBoundaryMesh &bm)
Return the indices of a processorCyclicPolyPatchs.
Calculates and applies the force necessary to maintain the specified mean velocity.
Namespace for OpenFOAM.