meshPhiCorrectInfoI.H
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) 2011-2022 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 "meshPhiCorrectInfo.H"
27 #include "fvMesh.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 inline void Foam::meshPhiCorrectInfo::validateFace() const
32 {
33  if (shape_ != shape::face)
34  {
36  << "Face data requested from a non-face correction info"
37  << exit(FatalError);
38  }
39 }
40 
41 
42 inline void Foam::meshPhiCorrectInfo::validateCell() const
43 {
44  if (shape_ != shape::cell)
45  {
47  << "Cell data requested from a non-cell correction info"
48  << exit(FatalError);
49  }
50 }
51 
52 
53 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
54 
56 :
57  shape_(shape::invalid),
58  delta_(NaN)
59 {}
60 
61 
63 :
64  shape_(s),
65  delta_(s == shape::invalid ? NaN : 0)
66 {}
67 
68 
69 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
70 
71 inline Foam::scalar Foam::meshPhiCorrectInfo::deltaPhi() const
72 {
73  validateFace();
74  return delta_;
75 }
76 
77 
78 inline Foam::scalar& Foam::meshPhiCorrectInfo::deltaPhi()
79 {
80  validateFace();
81  return delta_;
82 }
83 
84 
85 inline Foam::scalar Foam::meshPhiCorrectInfo::deltaDVdtError() const
86 {
87  validateCell();
88  return delta_;
89 }
90 
91 
93 {
94  validateCell();
95  return delta_;
96 }
97 
98 
99 template<class TrackingData>
100 inline bool Foam::meshPhiCorrectInfo::valid(TrackingData& td) const
101 {
102  return shape_ != shape::invalid;
103 }
104 
105 
106 template<class TrackingData>
108 (
109  const fvMesh&,
110  const meshPhiCorrectInfo& l,
111  const scalar tol,
112  TrackingData& td
113 ) const
114 {
115  return true;
116 }
117 
118 
119 template<class TrackingData>
121 (
122  const fvPatch& patch,
123  const label patchFacei,
124  const transformer& transform,
125  TrackingData& td
126 )
127 {}
128 
129 
130 template<class TrackingData>
132 (
133  const fvMesh& mesh,
134  const label thisCelli,
135  const labelPair& neighbourPatchAndFacei,
136  const meshPhiCorrectInfo& neighbourInfo,
137  const scalar tol,
138  TrackingData& td
139 )
140 {
141  const label neighbourPatchi = neighbourPatchAndFacei.first();
142  const label neighbourFacei = neighbourPatchAndFacei.second();
143 
144  const meshPhiPreCorrectInfo& thisPci =
145  td.cellPci_[thisCelli];
146  const meshPhiPreCorrectInfo& neighbourPci =
147  neighbourPatchi == -1
148  ? td.internalFacePci_[neighbourFacei]
149  : td.patchFacePci_[neighbourPatchi][neighbourFacei];
150 
151  if (!valid(td))
152  {
153  shape_ = shape::cell;
154  deltaDVdtError() = 0;
155  }
156 
157  if (thisPci.layer() < neighbourPci.layer())
158  {
159  const label s =
160  neighbourPatchi == -1
161  ? (mesh.owner()[neighbourFacei] == thisCelli ? +1 : -1)
162  : +1;
163 
164  deltaDVdtError() -= s*neighbourInfo.deltaPhi();
165 
166  return true;
167  }
168  else
169  {
170  return false;
171  }
172 }
173 
174 
175 template<class TrackingData>
177 (
178  const fvMesh& mesh,
179  const labelPair& thisPatchAndFacei,
180  const label neighbourCelli,
181  const meshPhiCorrectInfo& neighbourInfo,
182  const scalar tol,
183  TrackingData& td
184 )
185 {
186  const label thisPatchi = thisPatchAndFacei.first();
187  const label thisFacei = thisPatchAndFacei.second();
188 
189  const meshPhiPreCorrectInfo& thisPci =
190  thisPatchi == -1
191  ? td.internalFacePci_[thisFacei]
192  : td.patchFacePci_[thisPatchi][thisFacei];
193  const meshPhiPreCorrectInfo& neighbourPci =
194  td.cellPci_[neighbourCelli];
195 
196  if (!valid(td))
197  {
198  shape_ = shape::face;
199  deltaPhi() = 0;
200  }
201 
202  if (thisPci.layer() < neighbourPci.layer())
203  {
204  const label s =
205  thisPatchi == -1
206  ? (mesh.owner()[thisFacei] == neighbourCelli ? +1 : -1)
207  : +1;
208 
209  deltaPhi() =
210  s*thisPci.weight()/neighbourPci.weight()
211  *(td.dVdtError_[neighbourCelli] + neighbourInfo.deltaDVdtError());
212 
213  return true;
214  }
215  else
216  {
217  return false;
218  }
219 }
220 
221 
222 template<class TrackingData>
224 (
225  const fvMesh& mesh,
226  const labelPair& thisPatchAndFacei,
227  const meshPhiCorrectInfo& neighbourInfo,
228  const scalar tol,
229  TrackingData& td
230 )
231 {
232  if (!valid(td))
233  {
234  shape_ = shape::face;
235  deltaPhi() = 0;
236  }
237 
238  deltaPhi() = - neighbourInfo.deltaPhi();
239 
240  return true;
241 }
242 
243 
244 template<class TrackingData>
246 (
247  const meshPhiCorrectInfo& rhs,
248  TrackingData& td
249 ) const
250 {
251  return operator==(rhs);
252 }
253 
254 
255 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
256 
257 inline bool Foam::meshPhiCorrectInfo::operator==
258 (
259  const Foam::meshPhiCorrectInfo& rhs
260 ) const
261 {
262  return shape_ == rhs.shape_ && delta_ == rhs.delta_;
263 }
264 
265 
266 inline bool Foam::meshPhiCorrectInfo::operator!=
267 (
268  const Foam::meshPhiCorrectInfo& rhs
269 ) const
270 {
271  return !(*this == rhs);
272 }
273 
274 
275 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
276 
278 {
279  return os << static_cast<label>(l.shape_) << token::SPACE << l.delta_;
280 }
281 
282 
284 {
285  label s;
286  is >> s >> l.delta_;
287  l.shape_ = static_cast<meshPhiCorrectInfo::shape>(s);
288  return is;
289 }
290 
291 
292 // ************************************************************************* //
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
const Type & second() const
Return second.
Definition: Pair.H:110
const Type & first() const
Return first.
Definition: Pair.H:98
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:469
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
bool updateCell(const fvMesh &mesh, const label thisCelli, const labelPair &neighbourPatchAndFacei, const meshPhiCorrectInfo &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring face.
bool updateFace(const fvMesh &mesh, const labelPair &thisPatchAndFacei, const label neighbourCelli, const meshPhiCorrectInfo &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring cell.
scalar deltaPhi() const
Return the flux correction.
shape
Enumeration to define the mesh shape the info applies to.
meshPhiCorrectInfo()
Construct null.
bool equal(const meshPhiCorrectInfo &, TrackingData &td) const
Test equality.
bool sameGeometry(const fvMesh &mesh, const meshPhiCorrectInfo &, const scalar, TrackingData &td) const
Check for identical geometrical data. Used for checking.
scalar deltaDVdtError() const
Return the volume change rate error correction.
bool valid(TrackingData &td) const
Check whether the meshPhiCorrectInfo has been changed at all.
void transform(const fvPatch &patch, const label patchFacei, const transformer &transform, TrackingData &td)
Transform across an interface.
scalar weight() const
Return the weight.
label layer() const
Return the layer index.
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space.
Definition: transformer.H:84
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
bool valid(const PtrList< ModelType > &l)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Istream & operator>>(Istream &, pistonPointEdgeData &)
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
Ostream & operator<<(Ostream &os, const fvConstraints &constraints)
error FatalError