mappedFixedInternalValueFvPatchField.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-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 
27 #include "UIndirectList.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Type>
34 (
35  const fvPatch& p,
37 )
38 :
40 {}
41 
42 
43 template<class Type>
46 (
48  const fvPatch& p,
50  const fvPatchFieldMapper& mapper
51 )
52 :
53  mappedFixedValueFvPatchField<Type>(ptf, p, iF, mapper)
54 {}
55 
56 
57 template<class Type>
60 (
61  const fvPatch& p,
63  const dictionary& dict
64 )
65 :
67 {}
68 
69 
70 template<class Type>
73 (
75 )
76 :
78 {}
79 
80 
81 template<class Type>
84 (
87 )
88 :
90 {}
91 
92 
93 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
94 
95 template<class Type>
97 {
99 
100  if (this->updated())
101  {
102  return;
103  }
104 
105  // Since we're inside initEvaluate/evaluate there might be processor
106  // comms underway. Change the tag we use.
107  int oldTag = UPstream::msgType();
108  UPstream::msgType() = oldTag + 1;
109 
110  // Retrieve the neighbour values and assign to this patch boundary field
112 
113  // Get the coupling information from the mappedPatchBase
114  const mappedPatchBase& mpp =
115  refCast<const mappedPatchBase>(this->patch().patch());
116  const fvMesh& nbrMesh = refCast<const fvMesh>(mpp.sampleMesh());
117 
118  Field<Type> nbrIntFld;
119 
120  switch (mpp.mode())
121  {
122  case mappedPatchBase::NEARESTCELL:
123  {
125  (
126  "void mappedFixedValueFvPatchField<Type>::updateCoeffs()"
127  ) << "Cannot apply "
128  << mappedPatchBase::sampleModeNames_
129  [
130  mappedPatchBase::NEARESTCELL
131  ]
132  << " mapping mode for patch " << this->patch().name()
133  << exit(FatalError);
134 
135  break;
136  }
137  case mappedPatchBase::NEARESTPATCHFACE:
138  case mappedPatchBase::NEARESTPATCHFACEAMI:
139  {
140  const label samplePatchI = mpp.samplePolyPatch().index();
141  const fvPatchField<Type>& nbrPatchField =
142  this->sampleField().boundaryField()[samplePatchI];
143  nbrIntFld = nbrPatchField.patchInternalField();
144  mpp.distribute(nbrIntFld);
145 
146  break;
147  }
148  case mappedPatchBase::NEARESTFACE:
149  {
150  Field<Type> allValues(nbrMesh.nFaces(), pTraits<Type>::zero);
151 
152  const FieldType& nbrField = this->sampleField();
153 
154  forAll(nbrField.boundaryField(), patchI)
155  {
156  const fvPatchField<Type>& pf = nbrField.boundaryField()[patchI];
157  const Field<Type> pif(pf.patchInternalField());
158 
159  label faceStart = pf.patch().start();
160 
161  forAll(pf, faceI)
162  {
163  allValues[faceStart++] = pif[faceI];
164  }
165  }
166 
167  mpp.distribute(allValues);
168  nbrIntFld.transfer(allValues);
169 
170  break;
171  }
172  default:
173  {
174  FatalErrorIn("mappedFixedValueFvPatchField<Type>::updateCoeffs()")
175  << "Unknown sampling mode: " << mpp.mode()
176  << abort(FatalError);
177  }
178  }
179 
180  // Restore tag
181  UPstream::msgType() = oldTag;
182 
183  // Assign to (this) patch internal field its neighbour values
184  Field<Type>& intFld = const_cast<Field<Type>&>(this->internalField());
185  UIndirectList<Type>(intFld, this->patch().faceCells()) = nbrIntFld;
186 }
187 
188 
189 template<class Type>
191 (
192  Ostream& os
193 ) const
194 {
196 }
197 
198 
199 // ************************************************************************* //
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const polyMesh & sampleMesh() const
Get the region mesh.
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
label index() const
Return the index of this patch in the boundaryMesh.
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
Foam::fvPatchFieldMapper.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
runTime write()
dictionary dict
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const polyPatch & samplePolyPatch() const
Get the patch on the region.
This boundary condition maps the boundary and internal values of a neighbour patch field to the bound...
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:215
A List with indirect addressing.
Definition: fvMatrix.H:106
volScalarField & p
Definition: createFields.H:51
#define forAll(list, i)
Definition: UList.H:421
mappedFixedInternalValueFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Generic GeometricField class.
Traits class for primitives.
Definition: pTraits.H:50
const sampleMode & mode() const
What to sample.
error FatalError
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
This boundary condition maps the value at a set of cells or patch faces back to *this.
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:300
conserve internalField()+