36 template<
class Polynomial>
41 const scalar linearLimitFactor,
42 const scalar centralWeight
52 mesh, stencil, true, linearLimitFactor, centralWeight
54 coeffs_(
mesh.nFaces())
57 <<
"Constructing CentredFitSnGradData<Polynomial>" <<
nl;
61 DebugInfo <<
" Finished constructing polynomialFit data" <<
endl;
67 template<
class Polynomial>
73 const scalar deltaCoeff,
80 this->findFaceDirs(idir, jdir, kdir, facei);
84 wts[0] = this->centralWeight();
85 wts[1] = this->centralWeight();
116 Polynomial::addCoeffs(
B[ip], d, wts[ip], this->dim());
120 for (label i = 0; i <
B.m(); i++)
127 label stencilSize =
C.
size();
130 bool goodFit =
false;
131 for (
int iIt = 0; iIt < 8 && !goodFit; iIt++)
136 for (label i=0; i<stencilSize; i++)
138 coeffsi[i] = wts[1]*wts[i]*invB(1, i)/scale;
143 mag(wts[0]*wts[0]*invB(0, 0) - wLin)
144 < this->linearLimitFactor()*wLin)
145 && (
mag(wts[0]*wts[1]*invB(0, 1) - (1 - wLin)
146 ) < this->linearLimitFactor()*(1 - wLin))
147 && coeffsi[0] < 0 && coeffsi[1] > 0
148 &&
mag(coeffsi[0] + deltaCoeff) < 0.5*deltaCoeff
149 &&
mag(coeffsi[1] - deltaCoeff) < 0.5*deltaCoeff;
157 <<
"Cannot fit face " << facei <<
" iteration " << iIt
158 <<
" with sum of weights " <<
sum(coeffsi) <<
nl 159 <<
" Weights " << coeffsi <<
nl 160 <<
" Linear weights " << wLin <<
" " << 1 - wLin <<
nl 161 <<
" deltaCoeff " << deltaCoeff <<
nl 162 <<
" sing vals " << svd.S() <<
nl 163 <<
"Components of goodFit:\n" 164 <<
" wts[0]*wts[0]*invB(0, 0) = " 165 << wts[0]*wts[0]*invB(0, 0) <<
nl 166 <<
" wts[0]*wts[1]*invB(0, 1) = " 167 << wts[0]*wts[1]*invB(0, 1)
168 <<
" dim = " << this->dim() <<
endl;
173 for (label j = 0; j <
B.n(); j++)
179 for (label i = 0; i <
B.m(); i++)
190 coeffsi[0] += deltaCoeff;
191 coeffsi[1] -= deltaCoeff;
196 <<
"Could not fit face " << facei
197 <<
" Coefficients = " << coeffsi
198 <<
", reverting to uncorrected." <<
endl;
205 template<
class Polynomial>
213 List<List<point>> stencilPoints(
mesh.
nFaces());
214 this->stencil().collectData(
mesh.
C(), stencilPoints);
226 stencilPoints[facei],
243 label facei = pw.patch().start();
250 stencilPoints[facei],
fvsPatchField< scalar > fvsPatchScalarField
Graphite solid properties.
void size(const label n)
Older name for setAddressableSize.
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
GeometricBoundaryField< scalar, fvsPatchField, surfaceMesh > Boundary
Type of boundary fields.
label nFaces() const noexcept
Number of mesh faces.
CentredFitSnGradData(const fvMesh &mesh, const extendedCentredCellToFaceStencil &stencil, const scalar linearLimitFactor, const scalar centralWeight)
Construct from components.
#define forAll(list, i)
Loop across all elements in list.
virtual const surfaceScalarField & nonOrthDeltaCoeffs() const
Return reference to non-orthogonal cell-centre difference.
void setSize(const label n)
Alias for resize()
#define DebugInFunction
Report an information message using Foam::Info.
label nInternalFaces() const noexcept
Number of internal faces.
#define DebugInfo
Report an information message using Foam::Info.
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 calcFit()
Calculate the fit for all the faces.
const vectorField & faceCentres() const
#define WarningInFunction
Report a warning using Foam::Warning.
Polynomial templated on size (order):
Mesh data needed to do the Finite Volume discretisation.
Data for centred fit snGrad schemes.
const volVectorField & C() const
Return cell centres as volVectorField.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Data for the upwinded and centred polynomial fit interpolation schemes. The linearCorrection_ determi...
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
const volScalarField & p0
RectangularMatrix< scalar > scalarRectangularMatrix
static constexpr const zero Zero
Global zero (0)