findRefCell.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-2015 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 "findRefCell.H"
27 
28 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
29 
31 (
32  const volScalarField& field,
33  const volScalarField& fieldRef,
34  const dictionary& dict,
35  label& refCelli,
36  scalar& refValue,
37  const bool forceReference
38 )
39 {
40  if (fieldRef.needReference() || forceReference)
41  {
42  word refCellName = field.name() + "RefCell";
43  word refPointName = field.name() + "RefPoint";
44 
45  word refValueName = field.name() + "RefValue";
46 
47  if (dict.found(refCellName))
48  {
49  if (Pstream::master())
50  {
51  refCelli = readLabel(dict.lookup(refCellName));
52 
53  if (refCelli < 0 || refCelli >= field.mesh().nCells())
54  {
56  (
57  dict
58  ) << "Illegal master cellID " << refCelli
59  << ". Should be 0.." << field.mesh().nCells()
60  << exit(FatalIOError);
61  }
62  }
63  else
64  {
65  refCelli = -1;
66  }
67  }
68  else if (dict.found(refPointName))
69  {
70  point refPointi(dict.lookup(refPointName));
71 
72  // Try fast approximate search avoiding octree construction
73  refCelli = field.mesh().findCell(refPointi, polyMesh::FACE_PLANES);
74 
75  label hasRef = (refCelli >= 0 ? 1 : 0);
76  label sumHasRef = returnReduce<label>(hasRef, sumOp<label>());
77 
78  // If reference cell no found use octree search
79  // with cell tet-decompositoin
80  if (sumHasRef != 1)
81  {
82  refCelli = field.mesh().findCell(refPointi);
83 
84  hasRef = (refCelli >= 0 ? 1 : 0);
85  sumHasRef = returnReduce<label>(hasRef, sumOp<label>());
86  }
87 
88  if (sumHasRef != 1)
89  {
91  (
92  dict
93  ) << "Unable to set reference cell for field " << field.name()
94  << nl << " Reference point " << refPointName
95  << " " << refPointi
96  << " found on " << sumHasRef << " domains (should be one)"
97  << nl << exit(FatalIOError);
98  }
99  }
100  else
101  {
103  (
104  dict
105  ) << "Unable to set reference cell for field " << field.name()
106  << nl
107  << " Please supply either " << refCellName
108  << " or " << refPointName << nl << exit(FatalIOError);
109  }
110 
111  refValue = readScalar(dict.lookup(refValueName));
112  }
113 }
114 
115 
116 void Foam::setRefCell
117 (
118  const volScalarField& field,
119  const dictionary& dict,
120  label& refCelli,
121  scalar& refValue,
122  const bool forceReference
123 )
124 {
125  setRefCell(field, field, dict, refCelli, refValue, forceReference);
126 }
127 
128 
129 Foam::scalar Foam::getRefCellValue
130 (
131  const volScalarField& field,
132  const label refCelli
133 )
134 {
135  scalar refCellValue = (refCelli >= 0 ? field[refCelli] : 0.0);
136  return returnReduce(refCellValue, sumOp<scalar>());
137 }
138 
139 
140 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
bool needReference() const
Does the field need a reference level for solution.
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
scalar getRefCellValue(const volScalarField &field, const label refCelli)
Return the current value of field in the reference cell.
Definition: findRefCell.C:130
Find the reference cell nearest (in index) to the given cell but which is not on a cyclic...
A class for handling words, derived from string.
Definition: word.H:59
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
label readLabel(Istream &is)
Definition: label.H:64
static const char nl
Definition: Ostream.H:262
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
const Mesh & mesh() const
Return mesh.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
void setRefCell(const volScalarField &field, const volScalarField &fieldRef, const dictionary &dict, label &refCelli, scalar &refValue, const bool forceReference=false)
If the field fieldRef needs referencing find the reference cell nearest.
Definition: findRefCell.C:31
const word & name() const
Return name.
Definition: IOobject.H:260
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451
IOerror FatalIOError