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  "void Foam::setRefCell\n"
58  " (\n"
59  " const volScalarField&,\n"
60  " const volScalarField&,\n"
61  " const dictionary&,\n"
62  " label& scalar&,\n"
63  " bool\n"
64  ")",
65  dict
66  ) << "Illegal master cellID " << refCelli
67  << ". Should be 0.." << field.mesh().nCells()
68  << exit(FatalIOError);
69  }
70  }
71  else
72  {
73  refCelli = -1;
74  }
75  }
76  else if (dict.found(refPointName))
77  {
78  point refPointi(dict.lookup(refPointName));
79 
80  // Try fast approximate search avoiding octree construction
81  refCelli = field.mesh().findCell(refPointi, polyMesh::FACE_PLANES);
82 
83  label hasRef = (refCelli >= 0 ? 1 : 0);
84  label sumHasRef = returnReduce<label>(hasRef, sumOp<label>());
85 
86  // If reference cell no found use octree search
87  // with cell tet-decompositoin
88  if (sumHasRef != 1)
89  {
90  refCelli = field.mesh().findCell(refPointi);
91 
92  hasRef = (refCelli >= 0 ? 1 : 0);
93  sumHasRef = returnReduce<label>(hasRef, sumOp<label>());
94  }
95 
96  if (sumHasRef != 1)
97  {
99  (
100  "void Foam::setRefCell\n"
101  " (\n"
102  " const volScalarField&,\n"
103  " const volScalarField&,\n"
104  " const dictionary&,\n"
105  " label& scalar&,\n"
106  " bool\n"
107  " )",
108  dict
109  ) << "Unable to set reference cell for field " << field.name()
110  << nl << " Reference point " << refPointName
111  << " " << refPointi
112  << " found on " << sumHasRef << " domains (should be one)"
113  << nl << exit(FatalIOError);
114  }
115  }
116  else
117  {
119  (
120  "void Foam::setRefCell\n"
121  " (\n"
122  " const volScalarField&,\n"
123  " const volScalarField&,\n"
124  " const dictionary&,\n"
125  " label& scalar&,\n"
126  " bool\n"
127  " )",
128  dict
129  ) << "Unable to set reference cell for field " << field.name()
130  << nl
131  << " Please supply either " << refCellName
132  << " or " << refPointName << nl << exit(FatalIOError);
133  }
134 
135  refValue = readScalar(dict.lookup(refValueName));
136  }
137 }
138 
139 
140 void Foam::setRefCell
141 (
142  const volScalarField& field,
143  const dictionary& dict,
144  label& refCelli,
145  scalar& refValue,
146  const bool forceReference
147 )
148 {
149  setRefCell(field, field, dict, refCelli, refValue, forceReference);
150 }
151 
152 
153 Foam::scalar Foam::getRefCellValue
154 (
155  const volScalarField& field,
156  const label refCelli
157 )
158 {
159  scalar refCellValue = (refCelli >= 0 ? field[refCelli] : 0.0);
160  return returnReduce(refCellValue, sumOp<scalar>());
161 }
162 
163 
164 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
bool needReference() const
Does the field need a reference level for solution.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
#define readScalar
Definition: doubleScalar.C:38
A class for handling words, derived from string.
Definition: word.H:59
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
const Mesh & mesh() const
Return mesh.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
label readLabel(Istream &is)
Definition: label.H:64
static const char nl
Definition: Ostream.H:260
IOerror FatalIOError
Find the reference cell nearest (in index) to the given cell but which is not on a cyclic...
const word & name() const
Return name.
Definition: IOobject.H:260
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
scalar getRefCellValue(const volScalarField &field, const label refCelli)
Return the current value of field in the reference cell.
Definition: findRefCell.C:154
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325