Go to the documentation of this file.
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>
214 this->stencil().collectData(
mesh.C(), stencilPoints);
221 for (label facei = 0; facei <
mesh.nInternalFaces(); facei++)
226 stencilPoints[facei],
234 const surfaceScalarField::Boundary& bdC = dC.
boundaryField();
250 stencilPoints[facei],
const Cmpt & x() const
Access to the vector x component.
Data for the upwinded and centred polynomial fit interpolation schemes. The linearCorrection_ determi...
static constexpr const zero Zero
Global zero (0)
virtual 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...
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Cmpt & z() const
Access to the vector z component.
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.
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
scalarRectangularMatrix VSinvUt() const
Return the matrix product V S^(-1) U^T (the pseudo inverse)
#define DebugInFunction
Report an information message using Foam::Info.
Mesh data needed to do the Finite Volume discretisation.
Data for centred fit snGrad schemes.
const scalarDiagonalMatrix & S() const
Return the singular values.
const Cmpt & y() const
Access to the vector y component.
Polynomial templated on size (order):
#define DebugInfo
Report an information message using Foam::Info.
const fvPatch & patch() const
Return patch.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Singular value decomposition of a rectangular matrix.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const volScalarField & p0
Graphite solid properties.
void setSize(const label newSize)
Alias for resize(const label)
void calcFit()
Calculate the fit for all the faces.
virtual bool coupled() const
Return true if this patch field is coupled.
#define WarningInFunction
Report a warning using Foam::Warning.
const Boundary & boundaryField() const
Return const-reference to the boundary field.