acousticDampingSource.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-2021 OpenCFD Ltd.
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
29#include "fvMesh.H"
30#include "fvMatrices.H"
31#include "fvmSup.H"
35
36// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37
38namespace Foam
39{
40namespace fv
41{
44}
45}
46
47
48// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
49
51{
53
54 const vectorField& Cf = mesh_.C();
55
56 const scalar pi = constant::mathematical::pi;
57
58 forAll(cells_, i)
59 {
60 label celli = cells_[i];
61 scalar d = mag(Cf[celli] - x0_);
62
63 if (d < r1_)
64 {
65 blendFactor_[celli] = 0.0;
66 }
67 else if ((d >= r1_) && (d <= r2_))
68 {
69 blendFactor_[celli] = (1.0 - cos(pi*mag(d - r1_)/(r2_ - r1_)))/2.0;
70 }
71 }
72
74}
75
76
77// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
78
80(
81 const word& name,
82 const word& modelType,
83 const dictionary& dict,
84 const fvMesh& mesh
85)
86:
87 fv::cellSetOption(name, modelType, dict, mesh),
88 blendFactor_
89 (
91 (
93 (
94 name_ + ":blend",
95 mesh_.time().timeName(),
96 mesh_,
97 IOobject::NO_READ,
98 IOobject::NO_WRITE
99 ),
100 mesh_,
101 dimensionedScalar("blend0", dimless, 1.0),
103 )
104 ),
105 frequency_("frequency", dimless/dimTime, 0),
106 x0_(Zero),
107 r1_(0),
108 r2_(0),
109 URefName_("unknown-URef"),
110 w_(20)
111{
112 read(dict);
113}
114
115
116// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
117
119(
120 fvMatrix<vector>& eqn,
121 const label fieldI
122)
123{
124 const volVectorField& U = eqn.psi();
125 const volScalarField coeff(name_ + ":coeff", w_*frequency_*blendFactor_);
126 const auto& URef = mesh().lookupObject<volVectorField>(URefName_);
127
128 fvMatrix<vector> dampingEqn
129 (
130 fvm::Sp(coeff, U) - coeff*URef
131 );
132 eqn -= dampingEqn;
133}
134
135
137(
138 const volScalarField& rho,
139 fvMatrix<vector>& eqn,
140 const label fieldI
141)
142{
143 const volVectorField& U = eqn.psi();
144 const volScalarField coeff(name_ + ":coeff", w_*frequency_*blendFactor_);
145 const auto& URef = mesh().lookupObject<volVectorField>(URefName_);
146
147 fvMatrix<vector> dampingEqn
148 (
149 fvm::Sp(rho*coeff, U) - rho*coeff*URef
150 );
151 eqn -= dampingEqn;
152}
153
154
156(
157 const volScalarField& alpha,
158 const volScalarField& rho,
159 fvMatrix<vector>& eqn,
160 const label fieldI
161)
162{
163 const volVectorField& U = eqn.psi();
164 const volScalarField coeff(name_ + ":coeff", w_*frequency_*blendFactor_);
165 const auto& URef = mesh().lookupObject<volVectorField>(URefName_);
166
167 fvMatrix<vector> dampingEqn
168 (
169 fvm::Sp(alpha*rho*coeff, U) - alpha*rho*coeff*URef
170 );
171 eqn -= dampingEqn;
172}
173
174
176{
178 {
179 if (!coeffs_.readIfPresent("UNames", fieldNames_))
180 {
181 fieldNames_.resize(1);
182 fieldNames_.first() = coeffs_.getOrDefault<word>("U", "U");
183 }
184
186
187 coeffs_.readEntry("frequency", frequency_.value());
188 coeffs_.readEntry("URef", URefName_);
189 coeffs_.readCompat<vector>("origin", {{"centre", -1806}}, x0_);
190 coeffs_.readEntry("radius1", r1_);
191 coeffs_.readEntry("radius2", r2_);
192
193 if (coeffs_.readIfPresent("w", w_))
194 {
195 Info<< name_ << ": Setting stencil width to " << w_ << endl;
196 }
197
198 setBlendingFactor();
199
200 return true;
201 }
202
203 return false;
204}
205
206
207// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
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
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
const GeometricField< Type, fvPatchField, volMesh > & psi(const label i=0) const
Return psi.
Definition: fvMatrix.H:412
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const volVectorField & C() const
Return cell centres as volVectorField.
Applies sources on velocity (i.e. U) within a specified region to enable acoustic damping.
volScalarField blendFactor_
Blending factor [-].
virtual bool read(const dictionary &dict)
Read dictionary.
virtual void addSup(fvMatrix< vector > &eqn, const label fieldI)
Add implicit contribution to momentum equation.
void setBlendingFactor()
Helper function to set the blending factor.
Intermediate abstract class for handling cell-set options for the derived fvOptions.
labelList cells_
Set of cells to apply source to.
Base abstract class for handling finite volume options (i.e. fvOption).
Definition: fvOption.H:127
const fvMesh & mesh_
Reference to the mesh database.
Definition: fvOption.H:139
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
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
A special matrix type and solver, designed for finite volume solutions of scalar equations.
Calculate the finiteVolume matrix for implicit and explicit sources.
word timeName
Definition: getTimeIndex.H:3
constexpr scalar pi(M_PI)
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
Namespace for OpenFOAM.
const dimensionSet dimless
Dimensionless.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
dimensionedScalar cos(const dimensionedScalar &ds)
labelList fv(nPoints)
volScalarField & alpha
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333