33 template<
class Polynomial>
37 const extendedCentredCellToFaceStencil& stencil,
38 const scalar linearLimitFactor,
39 const scalar centralWeight
44 CentredFitSnGradData<Polynomial>,
45 extendedCentredCellToFaceStencil,
49 mesh, stencil,
true, linearLimitFactor, centralWeight
51 coeffs_(mesh.nFaces())
55 Info<<
"Contructing CentredFitSnGradData<Polynomial>" <<
endl;
62 Info<<
"CentredFitSnGradData<Polynomial>::CentredFitSnGradData() :" 63 <<
"Finished constructing polynomialFit data" 71 template<
class Polynomial>
77 const scalar deltaCoeff,
84 this->findFaceDirs(idir, jdir, kdir, facei);
88 wts[0] = this->centralWeight();
89 wts[1] = this->centralWeight();
92 point p0 = this->
mesh().faceCentres()[facei];
106 const vector p0p = p - p0;
120 Polynomial::addCoeffs(B[ip], d, wts[ip], this->dim());
124 for (
label i = 0; i < B.n(); i++)
134 bool goodFit =
false;
135 for (
int iIt = 0; iIt < 8 && !goodFit; iIt++)
139 for (
label i=0; i<stencilSize; i++)
141 coeffsi[i] = wts[1]*wts[i]*svd.
VSinvUt()[1][i]/scale;
147 < this->linearLimitFactor()*wLin)
148 && (
mag(wts[0]*wts[1]*svd.
VSinvUt()[0][1] - (1 - wLin)
149 ) < this->linearLimitFactor()*(1 - wLin))
150 && coeffsi[0] < 0 && coeffsi[1] > 0
151 &&
mag(coeffsi[0] + deltaCoeff) < 0.5*deltaCoeff
152 &&
mag(coeffsi[1] - deltaCoeff) < 0.5*deltaCoeff;
161 "CentredFitSnGradData<Polynomial>::calcFit" 162 "(const List<point>& C, const label facei" 163 ) <<
"Cannot fit face " << facei <<
" iteration " << iIt
164 <<
" with sum of weights " <<
sum(coeffsi) <<
nl 165 <<
" Weights " << coeffsi <<
nl 166 <<
" Linear weights " << wLin <<
" " << 1 - wLin <<
nl 167 <<
" deltaCoeff " << deltaCoeff <<
nl 168 <<
" sing vals " << svd.
S() <<
nl 169 <<
"Components of goodFit:\n" 170 <<
" wts[0]*wts[0]*svd.VSinvUt()[0][0] = " 171 << wts[0]*wts[0]*svd.
VSinvUt()[0][0] <<
nl 172 <<
" wts[0]*wts[1]*svd.VSinvUt()[0][1] = " 173 << wts[0]*wts[1]*svd.
VSinvUt()[0][1]
174 <<
" dim = " << this->dim() <<
endl;
179 for (
label j = 0; j < B.m(); j++)
185 for (
label i = 0; i < B.n(); i++)
196 coeffsi[0] += deltaCoeff;
197 coeffsi[1] -= deltaCoeff;
203 "CentredFitSnGradData<Polynomial>::calcFit(..)" 204 ) <<
"Could not fit face " << facei
205 <<
" Coefficients = " << coeffsi
206 <<
", reverting to uncorrected." <<
endl;
213 template<
class Polynomial>
222 this->stencil().collectData(mesh.
C(), stencilPoints);
234 stencilPoints[facei],
241 const surfaceScalarField::GeometricBoundaryField& bw = w.
boundaryField();
242 const surfaceScalarField::GeometricBoundaryField& bdC = dC.
boundaryField();
258 stencilPoints[facei],
Data for the upwinded and centred polynomial fit interpolation schemes. The linearCorrection_ determi...
virtual bool coupled() const
Return true if this patch field is coupled.
Polynomial templated on size (order):
Mesh data needed to do the Finite Volume discretisation.
dimensioned< scalar > mag(const dimensioned< Type > &)
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
const fvPatch & patch() const
Return patch.
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void size(const label)
Override size to be inconsistent with allocated storage.
void setSize(const label)
Reset size of List.
Ostream & endl(Ostream &os)
Add newline and flush stream.
#define WarningIn(functionName)
Report a warning using Foam::Warning.
const volVectorField & C() const
Return cell centres as volVectorField.
Singular value decomposition of a rectangular matrix.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const surfaceScalarField & nonOrthDeltaCoeffs() const
Return reference to non-orthogonal cell-centre difference.
label nInternalFaces() const
const scalarDiagonalMatrix & S() const
Return the singular values.
CentredFitSnGradData(const fvMesh &mesh, const extendedCentredCellToFaceStencil &stencil, const scalar linearLimitFactor, const scalar centralWeight)
Construct from components.
const scalarRectangularMatrix & VSinvUt() const
Return VSinvUt (the pseudo inverse)
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
label start() const
Return start label of this patch in the polyMesh face list.
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...