UpwindFitData.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 "UpwindFitData.H"
27 #include "surfaceFields.H"
28 #include "volFields.H"
29 #include "SVD.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
34 template<class Polynomial>
36 (
37  const fvMesh& mesh,
38  const extendedUpwindCellToFaceStencil& stencil,
39  const bool linearCorrection,
40  const scalar linearLimitFactor,
41  const scalar centralWeight
42 )
43 :
44  FitData
45  <
49  >
50  (
51  mesh, stencil, linearCorrection, linearLimitFactor, centralWeight
52  ),
53  owncoeffs_(mesh.nFaces()),
54  neicoeffs_(mesh.nFaces())
55 {
56  if (debug)
57  {
58  InfoInFunction << "Constructing UpwindFitData<Polynomial>" << endl;
59  }
60 
61  calcFit();
62 
63  if (debug)
64  {
65  Info<< " Finished constructing polynomialFit data" << endl;
66  }
67 }
68 
69 
70 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
71 
72 template<class Polynomial>
74 {
75  const fvMesh& mesh = this->mesh();
76 
77  const surfaceScalarField& w = mesh.surfaceInterpolation::weights();
79 
80  // Owner stencil weights
81  // ~~~~~~~~~~~~~~~~~~~~~
82 
83  // Get the cell/face centres in stencil order.
84  List<List<point>> stencilPoints(mesh.nFaces());
85  this->stencil().collectData
86  (
87  this->stencil().ownMap(),
88  this->stencil().ownStencil(),
89  mesh.C(),
90  stencilPoints
91  );
92 
93  // find the fit coefficients for every owner
94 
95  // Pout<< "-- Owner --" << endl;
96  for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
97  {
98  FitData
99  <
102  Polynomial
103  >::calcFit(owncoeffs_[facei], stencilPoints[facei], w[facei], facei);
104 
105  // Pout<< " facei:" << facei
106  // << " at:" << mesh.faceCentres()[facei] << endl;
107  // forAll(owncoeffs_[facei], i)
108  //{
109  // Pout<< " point:" << stencilPoints[facei][i]
110  // << "\tweight:" << owncoeffs_[facei][i]
111  // << endl;
112  //}
113  }
114 
115  forAll(bw, patchi)
116  {
117  const fvsPatchScalarField& pw = bw[patchi];
118 
119  if (pw.coupled())
120  {
121  label facei = pw.patch().start();
122 
123  forAll(pw, i)
124  {
125  FitData
126  <
127  UpwindFitData<Polynomial>,
128  extendedUpwindCellToFaceStencil,
129  Polynomial
130  >::calcFit
131  (
132  owncoeffs_[facei], stencilPoints[facei], pw[i], facei
133  );
134  facei++;
135  }
136  }
137  }
138 
139 
140  // Neighbour stencil weights
141  // ~~~~~~~~~~~~~~~~~~~~~~~~~
142 
143  // Note:reuse stencilPoints since is major storage
144  this->stencil().collectData
145  (
146  this->stencil().neiMap(),
147  this->stencil().neiStencil(),
148  mesh.C(),
149  stencilPoints
150  );
151 
152  // find the fit coefficients for every neighbour
153 
154  // Pout<< "-- Neighbour --" << endl;
155  for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
156  {
157  FitData
158  <
159  UpwindFitData<Polynomial>,
160  extendedUpwindCellToFaceStencil,
161  Polynomial
162  >::calcFit(neicoeffs_[facei], stencilPoints[facei], w[facei], facei);
163 
164  // Pout<< " facei:" << facei
165  // << " at:" << mesh.faceCentres()[facei] << endl;
166  // forAll(neicoeffs_[facei], i)
167  //{
168  // Pout<< " point:" << stencilPoints[facei][i]
169  // << "\tweight:" << neicoeffs_[facei][i]
170  // << endl;
171  //}
172  }
173 
174  forAll(bw, patchi)
175  {
176  const fvsPatchScalarField& pw = bw[patchi];
177 
178  if (pw.coupled())
179  {
180  label facei = pw.patch().start();
181 
182  forAll(pw, i)
183  {
184  FitData
185  <
186  UpwindFitData<Polynomial>,
187  extendedUpwindCellToFaceStencil,
188  Polynomial
189  >::calcFit
190  (
191  neicoeffs_[facei], stencilPoints[facei], pw[i], facei
192  );
193  facei++;
194  }
195  }
196  }
197 }
198 
199 
200 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Data for the upwinded and centred polynomial fit interpolation schemes. The linearCorrection_ determi...
Definition: FitData.H:57
Generic GeometricBoundaryField class.
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
Polynomial templated on size (order):
Definition: Polynomial.H:84
Data for the quadratic fit correction interpolation scheme to be used with upwind biased stencil.
Definition: UpwindFitData.H:62
UpwindFitData(const fvMesh &mesh, const extendedUpwindCellToFaceStencil &stencil, const bool linearCorrection, const scalar linearLimitFactor, const scalar centralWeight)
Construct from components.
Definition: UpwindFitData.C:36
Creates upwind stencil by shifting a centred stencil to upwind and downwind faces and optionally remo...
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const volVectorField & C() const
Return cell centres.
label nInternalFaces() const
label nFaces() const
label patchi
#define InfoInFunction
Report an information message using Foam::Info.
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
fvsPatchField< scalar > fvsPatchScalarField
Foam::surfaceFields.