inverseDistanceDiffusivity.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-2016 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 
28 #include "patchWave.H"
29 #include "HashSet.H"
30 #include "surfaceInterpolate.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(inverseDistanceDiffusivity, 0);
38 
40  (
41  motionDiffusivity,
42  inverseDistanceDiffusivity,
43  Istream
44  );
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
50 Foam::inverseDistanceDiffusivity::inverseDistanceDiffusivity
51 (
52  const fvMesh& mesh,
53  Istream& mdData
54 )
55 :
56  uniformDiffusivity(mesh, mdData),
57  patchNames_(mdData)
58 {
59  correct();
60 }
61 
62 
63 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
64 
66 {}
67 
68 
69 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
70 
71 Foam::tmp<Foam::scalarField> Foam::inverseDistanceDiffusivity::y() const
72 {
73  labelHashSet patchSet(mesh().boundaryMesh().patchSet(patchNames_));
74 
75  if (patchSet.size())
76  {
77  return tmp<scalarField>
78  (
79  new scalarField(patchWave(mesh(), patchSet, false).distance())
80  );
81  }
82  else
83  {
84  return tmp<scalarField>(new scalarField(mesh().nCells(), 1.0));
85  }
86 }
87 
88 
90 {
92  (
93  IOobject
94  (
95  "y",
96  mesh().time().timeName(),
97  mesh()
98  ),
99  mesh(),
100  dimless,
101  zeroGradientFvPatchScalarField::typeName
102  );
103  y_.primitiveFieldRef() = y();
105 
106  faceDiffusivity_ = 1.0/fvc::interpolate(y_);
107 }
108 
109 
110 // ************************************************************************* //
Uniform uniform finite volume mesh motion diffusivity.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
scalar distance(const vector &p1, const vector &p2)
Definition: curveTools.C:12
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
Macros for easy insertion into run-time selection tables.
scalar y
dynamicFvMesh & mesh
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Takes a set of patches to start MeshWave from. After construction holds distance at cells and distanc...
Definition: patchWave.H:56
word timeName
Definition: getTimeIndex.H:3
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil)-psi)*pSat, rhoMin);#1"/home/ubuntu/OpenFOAM-4.1/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:72
defineTypeNameAndDebug(combustionModel, 0)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual void correct()
Correct the motion diffusivity.
void correctBoundaryConditions()
Correct boundary field.
A class for managing temporary objects.
Definition: PtrList.H:54
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
Namespace for OpenFOAM.