porosityModel.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-2019 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 "porosityModel.H"
27 #include "volFields.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(porosityModel, 0);
34  defineRunTimeSelectionTable(porosityModel, mesh);
35 }
36 
37 
38 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
39 
41 {
42  scalar maxCmpt = cmptMax(resist.value());
43 
44  if (maxCmpt < 0)
45  {
47  << "Cannot have all resistances set to negative, resistance = "
48  << resist
49  << exit(FatalError);
50  }
51  else
52  {
53  maxCmpt = max(0, maxCmpt);
54  vector& val = resist.value();
55  for (label cmpt = 0; cmpt < vector::nComponents; cmpt++)
56  {
57  if (val[cmpt] < 0)
58  {
59  val[cmpt] *= -maxCmpt;
60  }
61  }
62  }
63 }
64 
65 
67 {
68  label index = 0;
69  if (!coordSys_.R().uniform())
70  {
71  index = i;
72  }
73  return index;
74 }
75 
76 
77 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
78 
79 Foam::porosityModel::porosityModel
80 (
81  const word& name,
82  const word& modelType,
83  const fvMesh& mesh,
84  const dictionary& dict,
85  const word& cellZoneName
86 )
87 :
89  (
90  IOobject
91  (
92  name,
93  mesh.time().timeName(),
94  mesh,
97  )
98  ),
99  name_(name),
100  mesh_(mesh),
101  dict_(dict),
102  coeffs_(dict.optionalSubDict(modelType + "Coeffs")),
103  active_(true),
104  zoneName_(cellZoneName),
105  cellZoneIDs_(),
106  coordSys_(coordinateSystem::New(mesh, coeffs_))
107 {
108  if (zoneName_ == word::null)
109  {
110  dict.readIfPresent("active", active_);
111  dict_.lookup("cellZone") >> zoneName_;
112  }
113 
114  cellZoneIDs_ = mesh_.cellZones().findIndices(zoneName_);
115 
116  Info<< " creating porous zone: " << zoneName_ << endl;
117 
118  bool foundZone = !cellZoneIDs_.empty();
119  reduce(foundZone, orOp<bool>());
120 
121  if (!foundZone && Pstream::master())
122  {
124  << "cannot find porous cellZone " << zoneName_
125  << exit(FatalError);
126  }
127 }
128 
129 
130 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
131 
133 {}
134 
135 
136 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
137 
139 {
140  if (!mesh_.upToDatePoints(*this))
141  {
142  calcTransformModelData();
143 
144  // set model up-to-date wrt points
145  mesh_.setUpToDatePoints(*this);
146  }
147 }
148 
149 
150 Foam::tmp<Foam::vectorField> Foam::porosityModel::porosityModel::force
151 (
152  const volVectorField& U,
153  const volScalarField& rho,
154  const volScalarField& mu
155 )
156 {
157  transformModelData();
158 
159  tmp<vectorField> tforce(new vectorField(U.size(), Zero));
160 
161  if (!cellZoneIDs_.empty())
162  {
163  this->calcForce(U, rho, mu, tforce.ref());
164  }
165 
166  return tforce;
167 }
168 
169 
171 {
172  if (cellZoneIDs_.empty())
173  {
174  return;
175  }
176 
177  transformModelData();
178  this->correct(UEqn);
179 }
180 
181 
183 (
184  fvVectorMatrix& UEqn,
185  const volScalarField& rho,
186  const volScalarField& mu
187 )
188 {
189  if (cellZoneIDs_.empty())
190  {
191  return;
192  }
193 
194  transformModelData();
195  this->correct(UEqn, rho, mu);
196 }
197 
198 
200 (
201  const fvVectorMatrix& UEqn,
202  volTensorField& AU,
203  bool correctAUprocBC
204 )
205 {
206  if (cellZoneIDs_.empty())
207  {
208  return;
209  }
210 
211  transformModelData();
212  this->correct(UEqn, AU);
213 
214  if (correctAUprocBC)
215  {
216  // Correct the boundary conditions of the tensorial diagonal to ensure
217  // processor boundaries are correctly handled when AU^-1 is interpolated
218  // for the pressure equation.
220  }
221 }
222 
223 
225 {
226  return true;
227 }
228 
229 
231 {
232  dict.readIfPresent("active", active_);
233 
234  coeffs_ = dict.optionalSubDict(type() + "Coeffs");
235 
236  dict.lookup("cellZone") >> zoneName_;
237  cellZoneIDs_ = mesh_.cellZones().findIndices(zoneName_);
238 
239  return true;
240 }
241 
242 
243 // ************************************************************************* //
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
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
virtual bool writeData(Ostream &os) const
Write.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual bool read()
Read object.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
label fieldIndex(const label index) const
Return label index.
Definition: porosityModel.C:66
void adjustNegativeResistance(dimensionedVector &resist)
Adjust negative resistance values to be multiplier of max value.
Definition: porosityModel.C:40
static const direction nComponents
Number of components in this vector space.
Definition: VectorSpace.H:99
dynamicFvMesh & mesh
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:766
A class for handling words, derived from string.
Definition: word.H:59
const Type & value() const
Return const reference to value.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
static const word null
An empty word.
Definition: word.H:77
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
static const zero Zero
Definition: zero.H:97
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual void transformModelData()
Transform the model data wrt mesh changes.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-7/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
defineTypeNameAndDebug(combustionModel, 0)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual ~porosityModel()
Destructor.
static autoPtr< coordinateSystem > New(const objectRegistry &obr, const dictionary &dict)
Select constructed from dictionary and objectRegistry.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
void correctBoundaryConditions()
Correct boundary field.
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:55
messageStream Info
Field< vector > vectorField
Specialisation of Field<T> for vector.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual void addResistance(fvVectorMatrix &UEqn)
Add resistance.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583