explicitPorositySource.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) 2012-2022 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 "explicitPorositySource.H"
27 #include "fvMesh.H"
28 #include "fvMatrices.H"
29 #include "porosityModel.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace fv
37 {
38  defineTypeNameAndDebug(explicitPorositySource, 0);
40  (
41  fvModel,
42  explicitPorositySource,
43  dictionary
44  );
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 void Foam::fv::explicitPorositySource::readCoeffs()
52 {
53  if (coeffs().found("UNames"))
54  {
55  UNames_ = wordList(coeffs().lookup("UNames"));
56  }
57  else
58  {
59  UNames_ = wordList(1, coeffs().lookupOrDefault<word>("U", "U"));
60  }
61 
62  porosityPtr_.reset
63  (
65  (
66  name(),
67  mesh(),
68  coeffs(),
69  set_.cellSetName()
70  ).ptr()
71  );
72 }
73 
74 
75 
76 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
77 
79 (
80  const word& name,
81  const word& modelType,
82  const dictionary& dict,
83  const fvMesh& mesh
84 )
85 :
86  fvModel(name, modelType, dict, mesh),
87  set_(coeffs(), mesh),
88  UNames_(),
89  porosityPtr_(nullptr)
90 {
91  readCoeffs();
92 }
93 
94 
95 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
96 
98 {
99  return UNames_;
100 }
101 
102 
104 (
105  fvMatrix<vector>& eqn,
106  const word& fieldName
107 ) const
108 {
109  fvMatrix<vector> porosityEqn(eqn.psi(), eqn.dimensions());
110  porosityPtr_->addResistance(porosityEqn);
111  eqn -= porosityEqn;
112 }
113 
114 
116 (
117  const volScalarField& rho,
118  fvMatrix<vector>& eqn,
119  const word& fieldName
120 ) const
121 {
122  fvMatrix<vector> porosityEqn(eqn.psi(), eqn.dimensions());
123  porosityPtr_->addResistance(porosityEqn);
124  eqn -= porosityEqn;
125 }
126 
127 
129 (
130  const volScalarField& alpha,
131  const volScalarField& rho,
132  fvMatrix<vector>& eqn,
133  const word& fieldName
134 ) const
135 {
136  fvMatrix<vector> porosityEqn(eqn.psi(), eqn.dimensions());
137  porosityPtr_->addResistance(porosityEqn);
138  eqn -= alpha*porosityEqn;
139 }
140 
141 
143 {
144  set_.movePoints();
145  return true;
146 }
147 
148 
150 {
151  set_.topoChange(map);
152 }
153 
154 
156 {
157  set_.mapMesh(map);
158 }
159 
160 
162 (
163  const polyDistributionMap& map
164 )
165 {
166  set_.distribute(map);
167 }
168 
169 
171 {
172  if (fvModel::read(dict))
173  {
174  set_.read(coeffs());
175  readCoeffs();
176  return true;
177  }
178  else
179  {
180  return false;
181  }
182 }
183 
184 
185 // ************************************************************************* //
defineTypeNameAndDebug(fixedTemperatureConstraint, 0)
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:290
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: fvModel.C:164
Finite volume model abstract base class.
Definition: fvModel.H:57
fvMesh & mesh
Macros for easy insertion into run-time selection tables.
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
virtual void addSup(fvMatrix< vector > &eqn, const word &fieldName) const
Add implicit contribution to momentum equation.
stressControl lookup("compactNormalStress") >> compactNormalStress
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
virtual bool movePoints()
Update for mesh motion.
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
virtual bool read(const dictionary &dict)
Read dictionary.
explicitPorositySource(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
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
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
addToRunTimeSelectionTable(fvConstraint, fixedTemperatureConstraint, dictionary)
List< word > wordList
A List of words.
Definition: fileName.H:54
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
A special matrix type and solver, designed for finite volume solutions of scalar equations.
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
static autoPtr< porosityModel > New(const word &name, const fvMesh &mesh, const dictionary &dict, const word &cellZoneName=word::null)
Selector.
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:50
bool found
Namespace for OpenFOAM.
const dimensionSet & dimensions() const
Definition: fvMatrix.H:295