All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
DarcyForchheimerTemplates.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 
26 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
27 
28 template<class RhoFieldType>
29 void Foam::porosityModels::DarcyForchheimer::apply
30 (
31  scalarField& Udiag,
32  vectorField& Usource,
33  const scalarField& V,
34  const RhoFieldType& rho,
35  const scalarField& mu,
36  const vectorField& U
37 ) const
38 {
39  forAll(cellZoneIDs_, zoneI)
40  {
41  const tensorField& dZones = D_[zoneI];
42  const tensorField& fZones = F_[zoneI];
43 
44  const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
45 
46  forAll(cells, i)
47  {
48  const label celli = cells[i];
49  const label j = this->fieldIndex(i);
50  const tensor Cd =
51  mu[celli]*dZones[j] + (rho[celli]*mag(U[celli]))*fZones[j];
52 
53  const scalar isoCd = tr(Cd);
54 
55  Udiag[celli] += V[celli]*isoCd;
56  Usource[celli] -= V[celli]*((Cd - I*isoCd) & U[celli]);
57  }
58  }
59 }
60 
61 
62 template<class RhoFieldType>
63 void Foam::porosityModels::DarcyForchheimer::apply
64 (
65  tensorField& AU,
66  const RhoFieldType& rho,
67  const scalarField& mu,
68  const vectorField& U
69 ) const
70 {
71  forAll(cellZoneIDs_, zoneI)
72  {
73  const tensorField& dZones = D_[zoneI];
74  const tensorField& fZones = F_[zoneI];
75 
76  const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zoneI]];
77 
78  forAll(cells, i)
79  {
80  const label celli = cells[i];
81  const label j = this->fieldIndex(i);
82  const tensor D = dZones[j];
83  const tensor F = fZones[j];
84 
85  AU[celli] += mu[celli]*D + (rho[celli]*mag(U[celli]))*F;
86  }
87  }
88 }
89 
90 
91 // ************************************************************************* //
const dimensionedScalar & F
Faraday constant: default SI units: [C/mol].
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
volVectorField vectorField(fieldObject, mesh)
label fieldIndex(const label index) const
Return label index.
Definition: porosityModel.C:66
const cellShapeList & cells
static const Identity< scalar > I
Definition: Identity.H:93
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:482
List< label > labelList
A List of labels.
Definition: labelList.H:56
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
volScalarField scalarField(fieldObject, mesh)
const fvMesh & mesh_
Reference to the mesh database.
Definition: porosityModel.H:76
labelList cellZoneIDs_
Cell zone IDs.
Definition: porosityModel.H:91
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
Tensor< scalar > tensor
Tensor of scalars.
Definition: tensor.H:51