cutPolyValueTemplates.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) 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 "cutPolyValue.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 template<class Type>
32 (
33  const edge& e,
34  const scalar lambda,
35  const Field<Type>& pPsis
36 )
37 {
38  return (1 - lambda)*pPsis[e[0]] + lambda*pPsis[e[1]];
39 }
40 
41 
42 template<class Type>
44 (
45  const edge& e,
46  const scalarField& pAlphas,
47  const scalar isoAlpha,
48  const Field<Type>& pPsis
49 )
50 {
51  return edgeCutValue(e, edgeCutLambda(e, pAlphas, isoAlpha), pPsis);
52 }
53 
54 
55 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
56 
57 template<class Type>
58 const Foam::Pair<Type>&
60 {
61  if (i != iAndCutPsis_.first())
62  {
63  const face& f = fValues_.f_;
64  const labelPair& fCut = fValues_.fCuts_[i];
65 
66  iAndCutPsis_.first() = i;
67 
68  forAll(iAndCutPsis_.second(), j)
69  {
70  iAndCutPsis_.second()[j] =
72  (
73  f.faceEdge(fCut[j]),
74  fValues_.pAlphas_,
75  fValues_.isoAlpha_,
76  fValues_.pPsis_
77  );
78  }
79  }
80 
81  return iAndCutPsis_.second();
82 }
83 
84 
85 template<class Type>
88 {
89  const face& f = fValues_.f_;
90  const List<labelPair>& fCuts = fValues_.fCuts_;
91  const labelPair& fCut0 = fCuts[i];
92  const labelPair& fCut1 = fCuts[(i + 1) % fCuts.size()];
93 
94  const label dfCut =
95  fValues_.below_ == separatedBelow
96  ? fCut0[1] - fCut0[0]
97  : fCut1[0] - fCut0[1];
98 
99  return 2 + ((dfCut + f.size()) % f.size());
100 }
101 
102 
103 template<class Type>
105 (
106  const label i,
107  const label j
108 ) const
109 {
110  const face& f = fValues_.f_;
111  const labelPair& fCut = fValues_.fCuts_[i];
112 
113  if (j < 2)
114  {
115  const label fCutj = fValues_.below_ == separatedBelow ? 1 - j : j;
116  return cutPsis(i)[fCutj];
117  }
118  else
119  {
120  const label fCutj = fValues_.below_ == separatedBelow ? 0 : 1;
121  return fValues_.pPsis_[f[(j - 2 + fCut[fCutj] + 1) % f.size()]];
122  }
123 }
124 
125 
126 template<class Type>
128 (
129  const label i,
130  const label j
131 )
132 {
133  const label iStar = i % cValues_.cCuts_.size();
134  const label jStar = j % cValues_.cCuts_[iStar].size();
135 
136  const label cei = cValues_.cCuts_[iStar][jStar];
137  const label cfi = cValues_.cAddr_.ceiToCfiAndFei()[cei][0][0];
138  const label fei = cValues_.cAddr_.ceiToCfiAndFei()[cei][0][1];
139 
140  return
142  (
143  cValues_.fs_[cValues_.c_[cfi]].faceEdge(fei),
144  cValues_.pAlphas_,
145  cValues_.isoAlpha_,
146  cValues_.pPsis_
147  );
148 }
149 
150 
151 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
152 
153 template<class Type>
155 (
156  const FaceCutValues<Type>& fValues,
157  const label i,
158  const label j
159 )
160 :
161  fValues_(fValues),
162  i_(i),
163  j_(j),
164  iAndCutPsis_(-1, Pair<Type>(Zero, Zero))
165 {}
166 
167 
168 template<class Type>
170 (
171  const face& f,
172  const List<labelPair>& fCuts,
173  const Field<Type>& pPsis,
174  const scalarField& pAlphas,
175  const scalar isoAlpha,
176  const bool below
177 )
178 :
179  f_(f),
180  fCuts_(fCuts),
181  pPsis_(pPsis),
182  pAlphas_(pAlphas),
183  isoAlpha_(isoAlpha),
184  below_(below)
185 {}
186 
187 
188 template<class Type>
190 (
191  const CellCutValues<Type>& cValues,
192  const label i,
193  const label j
194 )
195 :
196  cValues_(cValues),
197  i_(i),
198  j_(j),
199  psis_(psi(i, j), psi(i, j + 1))
200 {}
201 
202 
203 template<class Type>
205 (
206  const cell& c,
207  const cellEdgeAddressing& cAddr,
208  const labelListList& cCuts,
209  const faceList& fs,
210  const Field<Type>& pPsis,
211  const scalarField& pAlphas,
212  const scalar isoAlpha
213 )
214 :
215  c_(c),
216  cAddr_(cAddr),
217  cCuts_(cCuts),
218  fs_(fs),
219  pPsis_(pPsis),
220  pAlphas_(pAlphas),
221  isoAlpha_(isoAlpha)
222 {}
223 
224 
225 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
226 
227 template<class Type>
229 {
230  const label j = (j_ + 1) % size(i_);
231  const label i = i_ + (fValues_.below_ != separatedBelow && j == 0);
232  return i == fValues_.fCuts_.size() ? psi(0, 0) : psi(i, j);
233 }
234 
235 
236 template<class Type>
239 {
240  return const_iterator(*this, 0, 0);
241 }
242 
243 
244 template<class Type>
247 {
248  return const_iterator(*this, this->fCuts_.size(), 0);
249 }
250 
251 
252 template<class Type>
254 {
255  return psis_.second();
256 }
257 
258 
259 template<class Type>
262 {
263  return const_iterator(*this, 0, 0);
264 }
265 
266 
267 template<class Type>
270 {
271  return const_iterator(*this, this->cCuts_.size(), 0);
272 }
273 
274 
275 // * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
276 
277 template<class Type>
279 (
280  const const_iterator& it
281 ) const
282 {
283  return it.i_ == i_ && it.j_ == j_;
284 }
285 
286 
287 template<class Type>
289 (
290  const const_iterator& it
291 ) const
292 {
293  return !(it == *this);
294 }
295 
296 
297 template<class Type>
299 {
300  return psi(i_, j_);
301 }
302 
303 
304 template<class Type>
307 {
308  j_ = (j_ + 1) % size(i_);
309  i_ += j_ == 0;
310  return *this;
311 }
312 
313 
314 template<class Type>
317 {
318  const const_iterator it(*this);
319  ++ *this;
320  return it;
321 }
322 
323 
324 template<class Type>
326 (
327  const const_iterator& it
328 ) const
329 {
330  return it.i_ == i_ && it.j_ == j_;
331 }
332 
333 template<class Type>
335 (
336  const const_iterator& it
337 ) const
338 {
339  return !(it == *this);
340 }
341 
342 
343 template<class Type>
344 const Type&
346 {
347  return psis_.first();
348 }
349 
350 
351 template<class Type>
354 {
355  ++ j_;
356 
357  if (j_ == cValues_.cCuts_[i_].size())
358  {
359  j_ = 0;
360  ++ i_;
361  psis_.first() = psi(i_, j_);
362  }
363  else
364  {
365  psis_.first() = psis_.second();
366  }
367 
368  psis_.second() = psi(i_, j_ + 1);
369 
370  return *this;
371 }
372 
373 
374 template<class Type>
377 {
378  const const_iterator it(*this);
379  ++ *this;
380  return it;
381 }
382 
383 
384 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Pre-declare SubField and related Field type.
Definition: Field.H:83
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: Pair.H:65
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:60
const Type & next() const
Get the next value around the loop.
const_iterator(const CellCutValues< Type > &cValues, const label i, const label j)
Construct from components.
const_iterator begin() const
Get the beginning of the iteration over the values.
const_iterator end() const
Get the end of the iteration over the values.
CellCutValues(const cell &c, const cellEdgeAddressing &cAddr, const labelListList &cCuts, const faceList &fs, const Field< Type > &pPsis, const scalarField &pAlphas, const scalar isoAlpha)
Construct from components.
const_iterator(const FaceCutValues< Type > &fValues, const label i, const label j)
Construct from components.
Type next() const
Get the next value around the sub-face.
const_iterator end() const
Get the end of the iteration over the values.
const_iterator begin() const
Get the beginning of the iteration over the values.
FaceCutValues(const face &f, const List< labelPair > &fCuts, const Field< Type > &pPsis, const scalarField &pAlphas, const scalar isoAlpha, const bool below)
Construct from components.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:61
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:76
const volScalarField & psi
dimensionedScalar lambda(viscosity->lookup("lambda"))
const dimensionedScalar e
Elementary charge.
const dimensionedScalar c
Speed of light in a vacuum.
scalar edgeCutLambda(const edge &e, const scalarField &pAlphas, const scalar isoAlpha)
Get the local coordinate within an edge, given end point values and an.
Definition: cutPolyValueI.H:31
Type edgeCutValue(const edge &e, const scalar lambda, const Field< Type > &pPsis)
Linearly interpolate a value from the end points to the cut point of an.
static const bool separatedBelow
This flag determines which side of the iso-surface is separated into.
Definition: cutPoly.H:65
static const zero Zero
Definition: zero.H:97
Pair< label > labelPair
Label pair.
Definition: labelPair.H:48
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
labelList f(nPoints)