isotropicDamping.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) 2019-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 "isotropicDamping.H"
27 #include "fvMatrix.H"
28 #include "fvmSup.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace fv
36 {
37  defineTypeNameAndDebug(isotropicDamping, 0);
38  addToRunTimeSelectionTable(fvModel, isotropicDamping, dictionary);
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void Foam::fv::isotropicDamping::readCoeffs()
46 {
47  value_ =
49  (
50  value_.name(),
51  value_.dimensions(),
52  coeffs().lookup(value_.name())
53  );
54 }
55 
56 
57 void Foam::fv::isotropicDamping::add
58 (
59  const volScalarField::Internal& forceCoeff,
60  fvMatrix<vector>& eqn
61 ) const
62 {
63  eqn -= fvm::Sp(forceCoeff, eqn.psi());
64  eqn += forceCoeff*value_;
65 }
66 
67 
68 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
69 
71 (
72  const word& name,
73  const word& modelType,
74  const dictionary& dict,
75  const fvMesh& mesh
76 )
77 :
78  damping(name, modelType, dict, mesh),
79  value_("value", dimVelocity, vector::uniform(NaN))
80 {
81  readCoeffs();
82 }
83 
84 
85 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
86 
88 (
89  fvMatrix<vector>& eqn,
90  const word& fieldName
91 ) const
92 {
93  add(this->forceCoeff(), eqn);
94 }
95 
96 
98 (
99  const volScalarField& rho,
100  fvMatrix<vector>& eqn,
101  const word& fieldName
102 ) const
103 {
104  add(rho*forceCoeff(), eqn);
105 }
106 
107 
109 (
110  const volScalarField& alpha,
111  const volScalarField& rho,
112  fvMatrix<vector>& eqn,
113  const word& fieldName
114 ) const
115 {
116  add(alpha()*rho()*this->forceCoeff(), eqn);
117 }
118 
119 
121 {
122  if (damping::read(dict))
123  {
124  readCoeffs();
125  return true;
126  }
127  else
128  {
129  return false;
130  }
131 }
132 
133 
134 // ************************************************************************* //
virtual void addSup(fvMatrix< vector > &eqn, const word &fieldName) const
Source term to momentum equation.
defineTypeNameAndDebug(fixedTemperatureConstraint, 0)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
Macros for easy insertion into run-time selection tables.
virtual bool read(const dictionary &dict)
Read dictionary.
Definition: damping.C:201
isotropicDamping(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
virtual bool read(const dictionary &dict)
Read dictionary.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
void add(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
const dimensionSet dimVelocity
addToRunTimeSelectionTable(fvConstraint, fixedTemperatureConstraint, dictionary)
Base fvModel for damping functions.
Definition: damping.H:56
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
static Vector< scalar > uniform(const scalar &s)
Return a VectorSpace with all elements = s.
Definition: VectorSpaceI.H:168
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.