33 template<
class Form,
class ExtendedStencil,
class Polynomial>
37 const ExtendedStencil& stencil,
38 const bool linearCorrection,
39 const scalar linearLimitFactor,
40 const scalar centralWeight
45 linearCorrection_(linearCorrection),
46 linearLimitFactor_(linearLimitFactor),
47 centralWeight_(centralWeight),
48 #ifdef SPHERICAL_GEOMETRY 53 minSize_(Polynomial::nTerms(dim_))
56 if (linearLimitFactor <= small || linearLimitFactor > 3)
59 <<
"linearLimitFactor requested = " << linearLimitFactor
60 <<
" should be between zero and 3" 68 template<
class FitDataType,
class ExtendedStencil,
class Polynomial>
82 #ifndef SPHERICAL_GEOMETRY 111 kdir -= (idir & kdir)*idir;
113 scalar magk =
mag(kdir);
130 template<
class FitDataType,
class ExtendedStencil,
class Polynomial>
142 findFaceDirs(idir, jdir, kdir, facei);
146 wts[0] = centralWeight_;
147 if (linearCorrection_)
149 wts[1] = centralWeight_;
153 point p0 = this->
mesh().faceCentres()[facei];
171 d.
x() = (p - p0)&idir;
172 d.
y() = (p - p0)&jdir;
173 #ifndef SPHERICAL_GEOMETRY 174 d.
z() = (p - p0)&kdir;
187 Polynomial::addCoeffs
197 for (
label i = 0; i < B.m(); i++)
207 bool goodFit =
false;
208 for (
int iIt = 0; iIt < 8 && !goodFit; iIt++)
216 for (
label i=0; i<stencilSize; i++)
218 coeffsi[i] = wts[0]*wts[i]*invB(0, i);
219 if (
mag(coeffsi[i]) > maxCoeff)
221 maxCoeff =
mag(coeffsi[i]);
226 if (linearCorrection_)
229 (
mag(coeffsi[0] - wLin) < linearLimitFactor_*wLin)
230 && (
mag(coeffsi[1] - (1 - wLin)) < linearLimitFactor_*(1 - wLin))
237 (
mag(coeffsi[0] - 1.0) < linearLimitFactor_*1.0)
266 if (linearCorrection_)
271 for (
label j = 0; j < B.n(); j++)
277 for (
label i = 0; i < B.m(); i++)
287 if (linearCorrection_)
291 coeffsi[1] -= 1 - wLin;
304 <<
"Could not fit face " << facei
305 <<
" Weights = " << coeffsi
306 <<
", reverting to linear." <<
nl 307 <<
" Linear weights " << wLin <<
" " << 1 - wLin <<
endl;
315 template<
class FitDataType,
class ExtendedStencil,
class Polynomial>
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
virtual void calcFit()=0
Calculate the fit for all the faces.
errorManipArg< error, int > exit(error &err, const int errNo=1)
A face is a list of labels corresponding to mesh vertices.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Vector< scalar > vector
A scalar version of the templated Vector.
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
virtual const pointField & points() const
Return raw points.
const Vector< label > & geometricD() const
Return the vector of geometric directions in mesh.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Singular value decomposition of a rectangular matrix.
virtual const faceList & faces() const
Return raw faces.
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
scalarRectangularMatrix VSinvUt() const
Return the matrix product V S^(-1) U^T (the pseudo inverse)
const vectorField & faceCentres() const
void findFaceDirs(vector &idir, vector &jdir, vector &kdir, const label faci)
Find the normal direction (i) and j and k directions for face faci.
void setSize(const label)
Reset size of List.
#define WarningInFunction
Report a warning using Foam::Warning.
Mesh data needed to do the Finite Volume discretisation.
const vectorField & faceAreas() const
FitData(const fvMesh &mesh, const ExtendedStencil &stencil, const bool linearCorrection, const scalar linearLimitFactor, const scalar centralWeight)
Construct from components.
dimensioned< scalar > mag(const dimensioned< Type > &)
bool movePoints()
Recalculate weights (but not stencil) when the mesh moves.