inverseFaceDistanceDiffusivity.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2018 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
31 #include "HashSet.H"
32 #include "wallPoint.H"
33 #include "MeshWave.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(inverseFaceDistanceDiffusivity, 0);
40 
42  (
43  motionDiffusivity,
44  inverseFaceDistanceDiffusivity,
45  Istream
46  );
47 }
48 
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
52 Foam::inverseFaceDistanceDiffusivity::inverseFaceDistanceDiffusivity
53 (
54  const fvMesh& mesh,
55  Istream& mdData
56 )
57 :
58  uniformDiffusivity(mesh, mdData),
59  patchNames_(mdData)
60 {
61  correct();
62 }
63 
64 
65 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
66 
68 {}
69 
70 
71 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
72 
74 {
75  const polyBoundaryMesh& bdry = mesh().boundaryMesh();
76 
77  labelHashSet patchSet(bdry.size());
78 
79  label nPatchFaces = 0;
80 
81  for (const word& patchName : patchNames_)
82  {
83  const label patchi = bdry.findPatchID(patchName);
84 
85  if (patchi >= 0)
86  {
87  patchSet.insert(patchi);
88  nPatchFaces += bdry[patchi].size();
89  }
90  }
91 
92  List<wallPoint> faceDist(nPatchFaces);
93  labelList changedFaces(nPatchFaces);
94 
95  nPatchFaces = 0;
96 
97  for (const label patchi : patchSet)
98  {
99  const polyPatch& patch = bdry[patchi];
100 
101  const vectorField::subField fc(patch.faceCentres());
102 
103  forAll(fc, patchFacei)
104  {
105  changedFaces[nPatchFaces] = patch.start() + patchFacei;
106 
107  faceDist[nPatchFaces] = wallPoint(fc[patchFacei], 0);
108 
109  nPatchFaces++;
110  }
111  }
112  faceDist.setSize(nPatchFaces);
113  changedFaces.setSize(nPatchFaces);
114 
115  MeshWave<wallPoint> waveInfo
116  (
117  mesh(),
118  changedFaces,
119  faceDist,
120  mesh().globalData().nTotalCells()+1 // max iterations
121  );
122 
123  const List<wallPoint>& faceInfo = waveInfo.allFaceInfo();
124  const List<wallPoint>& cellInfo = waveInfo.allCellInfo();
125 
126  for (label facei=0; facei<mesh().nInternalFaces(); facei++)
127  {
128  scalar dist = faceInfo[facei].distSqr();
129 
130  faceDiffusivity_[facei] = 1.0/sqrt(dist);
131  }
132 
133  surfaceScalarField::Boundary& faceDiffusivityBf =
134  faceDiffusivity_.boundaryFieldRef();
135 
136  forAll(faceDiffusivityBf, patchi)
137  {
138  fvsPatchScalarField& bfld = faceDiffusivityBf[patchi];
139 
140  const labelUList& faceCells = bfld.patch().faceCells();
141 
142  if (patchSet.found(patchi))
143  {
144  forAll(bfld, i)
145  {
146  scalar dist = cellInfo[faceCells[i]].distSqr();
147  bfld[i] = 1.0/sqrt(dist);
148  }
149  }
150  else
151  {
152  const label start = bfld.patch().start();
153 
154  forAll(bfld, i)
155  {
156  scalar dist = faceInfo[start+i].distSqr();
157  bfld[i] = 1.0/sqrt(dist);
158  }
159  }
160  }
161 }
162 
163 
164 // ************************************************************************* //
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
Foam::wallPoint
Holds information regarding nearest wall point. Used in wall distance calculation.
Definition: wallPoint.H:65
Foam::MeshWave::allFaceInfo
const List< Type > & allFaceInfo() const
Get allFaceInfo.
Definition: MeshWave.H:125
wallPoint.H
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::inverseFaceDistanceDiffusivity::correct
virtual void correct()
Correct the motion diffusivity.
Definition: inverseFaceDistanceDiffusivity.C:73
Foam::polyBoundaryMesh
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Definition: polyBoundaryMesh.H:63
Foam::fvPatch::start
virtual label start() const
Return start label of this patch in the polyMesh face list.
Definition: fvPatch.H:173
Foam::uniformDiffusivity
Uniform uniform finite volume mesh motion diffusivity.
Definition: uniformDiffusivity.H:51
Foam::fvsPatchField
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:68
Foam::inverseFaceDistanceDiffusivity::~inverseFaceDistanceDiffusivity
virtual ~inverseFaceDistanceDiffusivity()
Destructor.
Definition: inverseFaceDistanceDiffusivity.C:67
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:444
Foam::HashSet< label, Hash< label > >
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::SubField
SubField is a Field obtained as a section of another Field, without its own allocation....
Definition: Field.H:64
Foam::Istream
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:61
MeshWave.H
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::MeshWave
FaceCellWave plus data.
Definition: MeshWave.H:59
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:222
correct
fvOptions correct(rho)
Foam::polyBoundaryMesh::findPatchID
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
Definition: polyBoundaryMesh.C:765
HashSet.H
inverseFaceDistanceDiffusivity.H
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::cellInfo
Holds information regarding type of cell. Used in inside/outside determination in cellClassification.
Definition: cellInfo.H:64
Foam::fvsPatchField::patch
const fvPatch & patch() const
Return patch.
Definition: fvsPatchField.H:281
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::primitiveMesh::nInternalFaces
label nInternalFaces() const noexcept
Number of internal faces.
Definition: primitiveMeshI.H:78
Foam::List
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: BitOps.H:63
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::UList< label >
Foam::fvPatch::faceCells
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:113
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::faceCells
Smooth ATC in cells next to a set of patches supplied by type.
Definition: faceCells.H:56
Foam::MeshWave::allCellInfo
const List< Type > & allCellInfo() const
Get allCellInfo.
Definition: MeshWave.H:131