limitVelocity.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) 2016-2017 OpenFOAM Foundation
9 Copyright (C) 2018-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "limitVelocity.H"
30#include "volFields.H"
32
33// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34
35namespace Foam
36{
37namespace fv
38{
41}
42}
43
44
45// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46
48(
49 const word& name,
50 const word& modelType,
51 const dictionary& dict,
52 const fvMesh& mesh
53)
54:
55 fv::cellSetOption(name, modelType, dict, mesh),
56 UName_(coeffs_.getOrDefault<word>("U", "U")),
57 max_(coeffs_.get<scalar>("max"))
58{
61}
62
63
64// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
65
67{
69 {
70 coeffs_.readEntry("max", max_);
71
72 return true;
73 }
74
75 return false;
76}
77
78
80{
81 const scalar maxSqrU = sqr(max_);
82
83 // Count nTotCells ourselves
84 // (maybe only applying on a subset)
85 label nCellsAbove(0);
86 const label nTotCells(returnReduce(cells_.size(), sumOp<label>()));
87
88 vectorField& Uif = U.primitiveFieldRef();
89
90 for (const label celli : cells_)
91 {
92 auto& Uval = Uif[celli];
93
94 const scalar magSqrUi = magSqr(Uval);
95
96 if (magSqrUi > maxSqrU)
97 {
98 Uval *= sqrt(maxSqrU/magSqrUi);
99 ++nCellsAbove;
100 }
101 }
102
103 // Handle boundaries in the case of 'all'
104
105 label nFacesAbove(0);
106 label nTotFaces(0);
107
109 {
110 for (fvPatchVectorField& Up : U.boundaryFieldRef())
111 {
112 if (!Up.fixesValue())
113 {
114 // Do not count patches that fix velocity themselves
115 nTotFaces += Up.size();
116
117 for (auto& Uval : Up)
118 {
119 const scalar magSqrUi = magSqr(Uval);
120
121 if (magSqrUi > maxSqrU)
122 {
123 Uval *= sqrt(maxSqrU/magSqrUi);
124 ++nFacesAbove;
125 }
126 }
127 }
128 }
129 }
130
131 // Percent, max 2 decimal places
132 const auto percent = [](scalar num, label denom) -> scalar
133 {
134 return (denom ? 1e-2*round(1e4*num/denom) : 0);
135 };
136
137
138 reduce(nCellsAbove, sumOp<label>());
139
140 // Report total numbers and percent
141 Info<< type() << ' ' << name_ << " Limited ";
142
143 Info<< nCellsAbove << " ("
144 << percent(nCellsAbove, nTotCells)
145 << "%) of cells";
146
147 reduce(nTotFaces, sumOp<label>());
148 reduce(nFacesAbove, sumOp<label>());
149 if (nTotFaces)
150 {
151 Info<< ", " << nFacesAbove << " ("
152 << percent(nFacesAbove, nTotFaces)
153 << "%) of faces";
154 }
155 Info<< ", with max limit " << max_ << endl;
156
157 if (nCellsAbove || nFacesAbove)
158 {
159 // We've changed internal values so give
160 // boundary conditions opportunity to correct
161 U.correctBoundaryConditions();
162 }
163}
164
165
166// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
virtual bool read()
Re-read model coefficients if they have changed.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Limits the maximum velocity magnitude to the specified max value.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Intermediate abstract class for handling cell-set options for the derived fvOptions.
bool useSubMesh() const noexcept
True if sub-selection should be used.
Corrects velocity field (i.e. T) within a specified region by applying a given maximum velocity magni...
word UName_
Name of operand velocity field.
Base abstract class for handling finite volume options (i.e. fvOption).
Definition: fvOption.H:127
wordList fieldNames_
Field names to apply source to - populated by derived models.
Definition: fvOption.H:148
void resetApplied()
Resize/reset applied flag list for all fieldNames_ entries.
Definition: fvOption.C:48
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
Namespace for OpenFOAM.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
messageStream Info
Information stream (stdout output on master, null elsewhere)
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensionedScalar sqrt(const dimensionedScalar &ds)
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
labelList fv(nPoints)
dictionary dict
volScalarField & e
Definition: createFields.H:11