FreeFOAM The Cross-Platform CFD Toolkit
volPointInterpolation.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) 1991-2010 OpenCFD Ltd.
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 "volPointInterpolation.H"
27 #include <finiteVolume/fvMesh.H>
28 #include <finiteVolume/volFields.H>
29 #include <OpenFOAM/pointFields.H>
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 defineTypeNameAndDebug(volPointInterpolation, 0);
42 
43 
44 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
45 
46 void volPointInterpolation::makeWeights()
47 {
48  if (debug)
49  {
50  Info<< "volPointInterpolation::makeWeights() : "
51  << "constructing weighting factors"
52  << endl;
53  }
54 
55  const pointField& points = mesh().points();
56  const labelListList& pointCells = mesh().pointCells();
57  const vectorField& cellCentres = mesh().cellCentres();
58 
59  // Allocate storage for weighting factors
60  pointWeights_.clear();
61  pointWeights_.setSize(points.size());
62 
63  forAll(pointWeights_, pointi)
64  {
65  pointWeights_[pointi].setSize(pointCells[pointi].size());
66  }
67 
68  pointScalarField sumWeights
69  (
70  IOobject
71  (
72  "volPointSumWeights",
74  mesh()
75  ),
77  dimensionedScalar("zero", dimless, 0)
78  );
79 
80  // Calculate inverse distances between cell centres and points
81  // and store in weighting factor array
82  forAll(points, pointi)
83  {
84  scalarList& pw = pointWeights_[pointi];
85  const labelList& pcp = pointCells[pointi];
86 
87  forAll(pcp, pointCelli)
88  {
89  pw[pointCelli] =
90  1.0/mag(points[pointi] - cellCentres[pcp[pointCelli]]);
91 
92  sumWeights[pointi] += pw[pointCelli];
93  }
94  }
95 
96  forAll(sumWeights.boundaryField(), patchi)
97  {
98  if (sumWeights.boundaryField()[patchi].coupled())
99  {
100  refCast<coupledPointPatchScalarField>
101  (sumWeights.boundaryField()[patchi]).initSwapAdd
102  (
103  sumWeights.internalField()
104  );
105  }
106  }
107 
108  forAll(sumWeights.boundaryField(), patchi)
109  {
110  if (sumWeights.boundaryField()[patchi].coupled())
111  {
112  refCast<coupledPointPatchScalarField>
113  (sumWeights.boundaryField()[patchi]).swapAdd
114  (
115  sumWeights.internalField()
116  );
117  }
118  }
119 
120  forAll(points, pointi)
121  {
122  scalarList& pw = pointWeights_[pointi];
123 
124  forAll(pw, pointCelli)
125  {
126  pw[pointCelli] /= sumWeights[pointi];
127  }
128  }
129 
130  if (debug)
131  {
132  Info<< "volPointInterpolation::makeWeights() : "
133  << "finished constructing weighting factors"
134  << endl;
135  }
136 }
137 
138 
139 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
140 
141 volPointInterpolation::volPointInterpolation(const fvMesh& vm)
142 :
144  boundaryInterpolator_(vm)
145 {
146  updateMesh();
147 }
148 
149 
150 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * //
151 
153 {}
154 
155 
156 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
157 
159 {
160  makeWeights();
161  boundaryInterpolator_.updateMesh();
162 }
163 
164 
166 {
167  makeWeights();
168  boundaryInterpolator_.movePoints();
169 
170  return true;
171 }
172 
173 
174 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
175 
176 } // End namespace Foam
177 
178 // ************************ vim: set sw=4 sts=4 et: ************************ //