singleCellFvMeshInterpolate.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) 2011-2018 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 "singleCellFvMesh.H"
29 #include "Time.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 
36 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
37 
38 template<class Type>
39 tmp<GeometricField<Type, fvPatchField, volMesh>> singleCellFvMesh::interpolate
40 (
42 ) const
43 {
44  // 1. Create the complete field with dummy patch fields
45  PtrList<fvPatchField<Type>> patchFields(vf.boundaryField().size());
46 
47  forAll(patchFields, patchi)
48  {
49  patchFields.set
50  (
51  patchi,
53  (
55  boundary()[patchi],
57  )
58  );
59  }
60 
61  // Create the complete field from the pieces
63  (
65  (
66  IOobject
67  (
68  vf.name(),
69  time().timeName(),
70  *this,
73  ),
74  *this,
75  vf.dimensions(),
76  Field<Type>(1, gAverage(vf)),
77  patchFields
78  )
79  );
81 
82 
83  // 2. Change the fvPatchFields to the correct type using a mapper
84  // constructor (with reference to the now correct internal field)
85 
87  Boundary& bf = resF.boundaryFieldRef();
88 
89  if (agglomerate())
90  {
92  {
93  const labelList& agglom = patchFaceAgglomeration_[patchi];
94  label nAgglom = max(agglom)+1;
95 
96  // Use inverse of agglomeration. This is from agglomeration to
97  // original (fine) mesh patch face.
98  labelListList coarseToFine(invertOneToMany(nAgglom, agglom));
99  inplaceReorder(patchFaceMap_[patchi], coarseToFine);
100  scalarListList coarseWeights(nAgglom);
101  forAll(coarseToFine, coarseI)
102  {
103  const labelList& fineFaces = coarseToFine[coarseI];
104  coarseWeights[coarseI] = scalarList
105  (
106  fineFaces.size(),
107  1.0/fineFaces.size()
108  );
109  }
110 
111  bf.set
112  (
113  patchi,
115  (
116  vf.boundaryField()[patchi],
117  boundary()[patchi],
118  resF(),
119  agglomPatchFieldMapper(coarseToFine, coarseWeights)
120  )
121  );
122  }
123  }
124  else
125  {
127  {
129 
130  bf.set
131  (
132  patchi,
134  (
135  vf.boundaryField()[patchi],
136  boundary()[patchi],
137  resF(),
139  )
140  );
141  }
142  }
143 
144  return tresF;
145 }
146 
147 
148 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
149 
150 } // End namespace Foam
151 
152 // ************************************************************************* //
void inplaceReorder(const labelUList &oldToNew, ListType &)
Inplace reorder the elements of a list.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const word & name() const
Return name.
Definition: IOobject.H:315
const Boundary & boundaryField() const
Return const-reference to the boundary field.
labelListList invertOneToMany(const label len, const labelUList &)
Invert one-to-many map. Unmapped elements will be size 0.
Definition: ListOps.C:67
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
Patch field mapper class for agglomerated meshes.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
tmp< GeometricField< Type, fvPatchField, volMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
Map volField. Internal field set to average, patch fields straight.
Generic GeometricField class.
const dimensionSet & dimensions() const
Return dimensions.
Pre-declare SubField and related Field type.
Definition: Field.H:56
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
This boundary condition is not designed to be evaluated; it is assumed that the value is assigned via...
faceListList boundary(nPatches)
Internal & ref()
Return a reference to the dimensioned internal field.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label patchi
Type gAverage(const FieldField< Field, Type > &f)
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
Namespace for OpenFOAM.