SemiImplicitSource.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 | Copyright (C) 2011-2015 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 "SemiImplicitSource.H"
27 #include "fvMesh.H"
28 #include "fvMatrices.H"
29 #include "DimensionedField.H"
30 #include "fvmSup.H"
31 
32 // * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * * //
33 
34 template<class Type>
36 (
37  IStringStream("(absolute specific)")()
38 );
39 
40 
41 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
42 
43 template<class Type>
46 (
47  const word& vmtName
48 ) const
49 {
50  forAll(volumeModeTypeNames_, i)
51  {
52  if (vmtName == volumeModeTypeNames_[i])
53  {
54  return volumeModeType(i);
55  }
56  }
57 
59  (
60  "SemiImplicitSource<Type>::volumeModeType"
61  "SemiImplicitSource<Type>::wordToVolumeModeType(const word&)"
62  ) << "Unknown volumeMode type " << vmtName
63  << ". Valid volumeMode types are:" << nl << volumeModeTypeNames_
64  << exit(FatalError);
65 
66  return volumeModeType(0);
67 }
68 
69 
70 template<class Type>
72 (
73  const volumeModeType& vmtType
74 ) const
75 {
76  if (vmtType > volumeModeTypeNames_.size())
77  {
78  return "UNKNOWN";
79  }
80  else
81  {
82  return volumeModeTypeNames_[vmtType];
83  }
84 }
85 
86 
87 template<class Type>
89 {
90  fieldNames_.setSize(dict.toc().size());
91  injectionRate_.setSize(fieldNames_.size());
92 
93  applied_.setSize(fieldNames_.size(), false);
94 
95  label i = 0;
96  forAllConstIter(dictionary, dict, iter)
97  {
98  fieldNames_[i] = iter().keyword();
99  dict.lookup(iter().keyword()) >> injectionRate_[i];
100  i++;
101  }
102 
103  // Set volume normalisation
104  if (volumeMode_ == vmAbsolute)
105  {
106  VDash_ = V_;
107  }
108 }
109 
110 
111 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
112 
113 template<class Type>
115 (
116  const word& name,
117  const word& modelType,
118  const dictionary& dict,
119  const fvMesh& mesh
120 )
121 :
122  cellSetOption(name, modelType, dict, mesh),
123  volumeMode_(vmAbsolute),
124  VDash_(1.0),
125  injectionRate_()
126 {
127  read(dict);
128 }
129 
130 
131 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
132 
133 template<class Type>
135 (
136  fvMatrix<Type>& eqn,
137  const label fieldI
138 )
139 {
140  if (debug)
141  {
142  Info<< "SemiImplicitSource<" << pTraits<Type>::typeName
143  << ">::addSup for source " << name_ << endl;
144  }
145 
147 
149  (
150  IOobject
151  (
152  name_ + fieldNames_[fieldI] + "Su",
153  mesh_.time().timeName(),
154  mesh_,
155  IOobject::NO_READ,
156  IOobject::NO_WRITE
157  ),
158  mesh_,
160  (
161  "zero",
162  eqn.dimensions()/dimVolume,
164  ),
165  false
166  );
167 
168  UIndirectList<Type>(Su, cells_) = injectionRate_[fieldI].first()/VDash_;
169 
171  (
172  IOobject
173  (
174  name_ + fieldNames_[fieldI] + "Sp",
175  mesh_.time().timeName(),
176  mesh_,
177  IOobject::NO_READ,
178  IOobject::NO_WRITE
179  ),
180  mesh_,
182  (
183  "zero",
184  Su.dimensions()/psi.dimensions(),
185  0.0
186  ),
187  false
188  );
189 
190  UIndirectList<scalar>(Sp, cells_) = injectionRate_[fieldI].second()/VDash_;
191 
192  eqn += Su + fvm::SuSp(Sp, psi);
193 }
194 
195 
196 template<class Type>
198 (
199  const volScalarField& rho,
200  fvMatrix<Type>& eqn,
201  const label fieldI
202 )
203 {
204  if (debug)
205  {
206  Info<< "SemiImplicitSource<" << pTraits<Type>::typeName
207  << ">::addSup for source " << name_ << endl;
208  }
209 
210  return this->addSup(eqn, fieldI);
211 }
212 
213 
214 // ************************************************************************* //
Cell-set options abtract base class. Provides a base set of controls, e.g.
Definition: cellSetOption.H:71
word volumeModeTypeToWord(const volumeModeType &vtType) const
Helper function to convert from a volumeModeType to a word.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dimensionSet & dimensions() const
Definition: fvMatrix.H:286
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Calculate the matrix for implicit and explicit sources.
A class for handling words, derived from string.
Definition: word.H:59
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
tmp< GeometricField< Type, fvPatchField, volMesh > > Sp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:67
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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:68
messageStream Info
SemiImplicitSource(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
static const wordList volumeModeTypeNames_
Word list of volume mode type names.
tmp< GeometricField< Type, fvPatchField, volMesh > > Su(const GeometricField< Type, fvPatchField, volMesh > &su, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:44
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:281
Generic dimensioned Type class.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
A List with indirect addressing.
Definition: fvMatrix.H:106
volumeModeType wordToVolumeModeType(const word &vtName) const
Helper function to convert from a word to a volumeModeType.
wordList toc() const
Return the table of contents.
Definition: dictionary.C:707
volumeModeType
Enumeration for volume types.
#define forAll(list, i)
Definition: UList.H:421
tmp< GeometricField< Type, fvPatchField, volMesh > > SuSp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:102
const volScalarField & psi
Definition: createFields.H:24
const dimensionSet & dimensions() const
Return dimensions.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Generic GeometricField class.
Traits class for primitives.
Definition: pTraits.H:50
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
virtual void addSup(fvMatrix< Type > &eqn, const label fieldI)
Add explicit contribution to equation.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
bool read(const char *, int32_t &)
Definition: int32IO.C:87
void setFieldData(const dictionary &dict)
Set the local field data.
A special matrix type and solver, designed for finite volume solutions of scalar equations.