invDistLeastSquaresVectors.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) 2013 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 "leastSquaresVectors.H"
27 #include "volFields.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(leastSquaresVectors, 0);
34 }
35 
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
40 :
42  pVectors_
43  (
44  IOobject
45  (
46  "LeastSquaresP",
47  mesh_.pointsInstance(),
48  mesh_,
49  IOobject::NO_READ,
50  IOobject::NO_WRITE,
51  false
52  ),
53  mesh_,
55  ),
56  nVectors_
57  (
58  IOobject
59  (
60  "LeastSquaresN",
61  mesh_.pointsInstance(),
62  mesh_,
63  IOobject::NO_READ,
64  IOobject::NO_WRITE,
65  false
66  ),
67  mesh_,
69  )
70 {
71  calcLeastSquaresVectors();
72 }
73 
74 
75 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
76 
78 {}
79 
80 
81 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
82 
83 void Foam::leastSquaresVectors::calcLeastSquaresVectors()
84 {
85  if (debug)
86  {
87  Info<< "leastSquaresVectors::calcLeastSquaresVectors() :"
88  << "Calculating least square gradient vectors"
89  << endl;
90  }
91 
92  const fvMesh& mesh = mesh_;
93 
94  // Set local references to mesh data
95  const labelUList& owner = mesh_.owner();
96  const labelUList& neighbour = mesh_.neighbour();
97 
98  const volVectorField& C = mesh.C();
99 
100  // Set up temporary storage for the dd tensor (before inversion)
102 
103  forAll(owner, facei)
104  {
105  label own = owner[facei];
106  label nei = neighbour[facei];
107 
108  vector d = C[nei] - C[own];
109  symmTensor wdd = sqr(d)/magSqr(d);
110  dd[own] += wdd;
111  dd[nei] += wdd;
112  }
113 
114 
115  surfaceVectorField::GeometricBoundaryField& blsP =
116  pVectors_.boundaryField();
117 
118  forAll(blsP, patchi)
119  {
120  const fvsPatchVectorField& patchLsP = blsP[patchi];
121 
122  const fvPatch& p = patchLsP.patch();
123  const labelUList& faceCells = p.patch().faceCells();
124 
125  // Build the d-vectors
126  vectorField pd(p.delta());
127 
128  forAll(pd, patchFacei)
129  {
130  const vector& d = pd[patchFacei];
131 
132  dd[faceCells[patchFacei]] += sqr(d)/magSqr(d);
133  }
134  }
135 
136 
137  // Invert the dd tensor
138  const symmTensorField invDd(inv(dd));
139 
140 
141  // Revisit all faces and calculate the pVectors_ and nVectors_ vectors
142  forAll(owner, facei)
143  {
144  label own = owner[facei];
145  label nei = neighbour[facei];
146 
147  vector d = C[nei] - C[own];
148 
149  pVectors_[facei] = (invDd[own] & d)/magSqr(d);
150  nVectors_[facei] = -(invDd[nei] & d)/magSqr(d);
151  }
152 
153  forAll(blsP, patchi)
154  {
155  fvsPatchVectorField& patchLsP = blsP[patchi];
156 
157  const fvPatch& p = patchLsP.patch();
158  const labelUList& faceCells = p.faceCells();
159 
160  // Build the d-vectors
161  vectorField pd(p.delta());
162 
163  forAll(pd, patchFacei)
164  {
165  const vector& d = pd[patchFacei];
166 
167  patchLsP[patchFacei] = (invDd[faceCells[patchFacei]] & d)/magSqr(d);
168  }
169  }
170 
171  if (debug)
172  {
173  Info<< "leastSquaresVectors::calcLeastSquaresVectors() :"
174  << "Finished calculating least square gradient vectors"
175  << endl;
176  }
177 }
178 
179 
181 {
182  calcLeastSquaresVectors();
183  return true;
184 }
185 
186 
187 // ************************************************************************* //
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
const fvPatch & patch() const
Return patch.
leastSquaresVectors(const fvMesh &)
Construct given an fvMesh.
Templated 3D symmetric tensor derived from VectorSpace adding construction from 6 components...
Definition: SymmTensor.H:53
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual tmp< vectorField > delta() const
Return cell-centre to face-centre vector.
Definition: fvPatch.C:142
virtual ~leastSquaresVectors()
Destructor.
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:316
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
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:282
label nCells() const
Graphite solid properties.
Definition: C.H:57
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:93
messageStream Info
dynamicFvMesh & mesh
static const SymmTensor zero
Definition: SymmTensor.H:77
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
Namespace for OpenFOAM.
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:143
Least-squares gradient scheme vectors.
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
const volVectorField & C() const
Return cell centres as volVectorField.
volScalarField & p
Definition: createFields.H:51
#define forAll(list, i)
Definition: UList.H:421
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:47
label patchi
virtual bool movePoints()
Delete the least square vectors when the mesh moves.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Definition: MeshObject.H:81
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:288
dimensionedSymmTensor sqr(const dimensionedVector &dv)
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
defineTypeNameAndDebug(combustionModel, 0)