powerLaw.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-2018 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 
27 #include "powerLaw.H"
28 #include "geometricOneField.H"
29 #include "fvMatrices.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  namespace porosityModels
36  {
37  defineTypeNameAndDebug(powerLaw, 0);
38  addToRunTimeSelectionTable(porosityModel, powerLaw, mesh);
39  }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
45 Foam::porosityModels::powerLaw::powerLaw
46 (
47  const word& name,
48  const word& modelType,
49  const fvMesh& mesh,
50  const dictionary& dict,
51  const word& cellZoneName
52 )
53 :
54  porosityModel(name, modelType, mesh, dict, cellZoneName),
55  C0_(readScalar(coeffs_.lookup("C0"))),
56  C1_(readScalar(coeffs_.lookup("C1"))),
57  rhoName_(coeffs_.lookupOrDefault<word>("rho", "rho"))
58 {}
59 
60 
61 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
62 
64 {}
65 
66 
67 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
68 
70 {
71  // nothing to be transformed
72 }
73 
74 
76 (
77  const volVectorField& U,
78  const volScalarField& rho,
79  const volScalarField& mu,
80  vectorField& force
81 ) const
82 {
83  scalarField Udiag(U.size(), 0.0);
84  const scalarField& V = mesh_.V();
85 
86  apply(Udiag, V, rho, U);
87 
88  force = Udiag*U;
89 }
90 
91 
93 (
94  fvVectorMatrix& UEqn
95 ) const
96 {
97  const volVectorField& U = UEqn.psi();
98  const scalarField& V = mesh_.V();
99  scalarField& Udiag = UEqn.diag();
100 
101  if (UEqn.dimensions() == dimForce)
102  {
103  const volScalarField& rho = mesh_.lookupObject<volScalarField>
104  (
105  IOobject::groupName(rhoName_, U.group())
106  );
107 
108  apply(Udiag, V, rho, U);
109  }
110  else
111  {
112  apply(Udiag, V, geometricOneField(), U);
113  }
114 }
115 
116 
118 (
119  fvVectorMatrix& UEqn,
120  const volScalarField& rho,
121  const volScalarField& mu
122 ) const
123 {
124  const vectorField& U = UEqn.psi();
125  const scalarField& V = mesh_.V();
126  scalarField& Udiag = UEqn.diag();
127 
128  apply(Udiag, V, rho, U);
129 }
130 
131 
133 (
134  const fvVectorMatrix& UEqn,
135  volTensorField& AU
136 ) const
137 {
138  const volVectorField& U = UEqn.psi();
139 
140  if (UEqn.dimensions() == dimForce)
141  {
142  const volScalarField& rho = mesh_.lookupObject<volScalarField>
143  (
144  IOobject::groupName(rhoName_, U.group())
145  );
146 
147  apply(AU, rho, U);
148  }
149  else
150  {
151  apply(AU, geometricOneField(), U);
152  }
153 }
154 
155 
157 {
158  os << indent << name_ << endl;
159  dict_.write(os);
160 
161  return true;
162 }
163 
164 
165 // ************************************************************************* //
static word group(const word &name)
Return group (extension part of name)
Definition: IOobject.C:176
defineTypeNameAndDebug(powerLawLopesdaCosta, 0)
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:226
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
addToRunTimeSelectionTable(porosityModel, powerLawLopesdaCosta, mesh)
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:284
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
Macros for easy insertion into run-time selection tables.
virtual void calcTransformModelData()
Transform the model data wrt mesh changes.
Definition: powerLaw.C:69
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
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
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const dimensionSet dimForce
virtual void calcForce(const volVectorField &U, const volScalarField &rho, const volScalarField &mu, vectorField &force) const
Calculate the porosity force.
Definition: powerLaw.C:76
virtual ~powerLaw()
Destructor.
Definition: powerLaw.C:63
virtual void correct(fvVectorMatrix &UEqn) const
Add resistance.
Definition: powerLaw.C:93
bool writeData(Ostream &os) const
Write.
Definition: powerLaw.C:156
U
Definition: pEqn.H:72
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A special matrix type and solver, designed for finite volume solutions of scalar equations.
virtual Ostream & write(const token &)=0
Write next token to stream.
scalarField & diag()
Definition: lduMatrix.C:186
Top level model for porosity models.
Definition: porosityModel.H:55
Namespace for OpenFOAM.
const dimensionSet & dimensions() const
Definition: fvMatrix.H:289