35 template<
class Form,
class ExtendedStencil,
class Polynomial>
39 const ExtendedStencil& stencil,
40 const bool linearCorrection,
41 const scalar linearLimitFactor,
42 const scalar centralWeight
47 linearCorrection_(linearCorrection),
48 linearLimitFactor_(linearLimitFactor),
49 centralWeight_(centralWeight),
50 #ifdef SPHERICAL_GEOMETRY
53 dim_(
mesh.nGeometricD()),
58 if (linearLimitFactor <= SMALL || linearLimitFactor > 3)
61 <<
"linearLimitFactor requested = " << linearLimitFactor
62 <<
" should be between zero and 3" 70 template<
class FitDataType,
class ExtendedStencil,
class Polynomial>
83 #ifndef SPHERICAL_GEOMETRY 84 if (
mesh.nGeometricD() <= 2)
86 if (
mesh.geometricD()[0] == -1)
90 else if (
mesh.geometricD()[1] == -1)
101 const face&
f =
mesh.faces()[facei];
102 kdir =
mesh.points()[
f[0]] -
mesh.faceCentres()[facei];
106 kdir =
mesh.faceCentres()[facei];
109 if (
mesh.nGeometricD() == 3)
112 kdir -= (idir & kdir)*idir;
114 scalar magk =
mag(kdir);
131 template<
class FitDataType,
class ExtendedStencil,
class Polynomial>
143 findFaceDirs(idir, jdir, kdir, facei);
147 wts[0] = centralWeight_;
148 if (linearCorrection_)
150 wts[1] = centralWeight_;
172 d.x() = (
p -
p0)&idir;
173 d.y() = (
p -
p0)&jdir;
174 #ifndef SPHERICAL_GEOMETRY 175 d.z() = (
p -
p0)&kdir;
188 Polynomial::addCoeffs
198 for (label i = 0; i <
B.m(); i++)
205 label stencilSize =
C.size();
208 bool goodFit =
false;
209 for (
int iIt = 0; iIt < 8 && !goodFit; iIt++)
217 for (label i=0; i<stencilSize; i++)
219 coeffsi[i] = wts[0]*wts[i]*invB(0, i);
220 if (
mag(coeffsi[i]) > maxCoeff)
222 maxCoeff =
mag(coeffsi[i]);
227 if (linearCorrection_)
230 (
mag(coeffsi[0] - wLin) < linearLimitFactor_*wLin)
231 && (
mag(coeffsi[1] - (1 - wLin)) < linearLimitFactor_*(1 - wLin))
238 (
mag(coeffsi[0] - 1.0) < linearLimitFactor_*1.0)
267 if (linearCorrection_)
272 for (label j = 0; j <
B.n(); j++)
278 for (label i = 0; i <
B.m(); i++)
288 if (linearCorrection_)
292 coeffsi[1] -= 1 - wLin;
305 <<
"Could not fit face " << facei
306 <<
" Weights = " << coeffsi
307 <<
", reverting to linear." <<
nl 308 <<
" Linear weights " << wLin <<
" " << 1 - wLin <<
endl;
316 template<
class FitDataType,
class ExtendedStencil,
class Polynomial>
Graphite solid properties.
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
virtual void calcFit()=0
Calculate the fit for all the faces.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
#define forAll(list, i)
Loop across all elements in list.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
void setSize(const label n)
Alias for resize()
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
void findFaceDirs(vector &idir, vector &jdir, vector &kdir, const label faci)
Find the normal direction (i) and j and k directions for face faci.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.
Polynomial templated on size (order):
Mesh data needed to do the Finite Volume discretisation.
FitData(const fvMesh &mesh, const ExtendedStencil &stencil, const bool linearCorrection, const scalar linearLimitFactor, const scalar centralWeight)
Construct from components.
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
const volScalarField & p0
bool movePoints()
Recalculate weights (but not stencil) when the mesh moves.
RectangularMatrix< scalar > scalarRectangularMatrix
static constexpr const zero Zero
Global zero (0)