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.
void setSize(const label n)
Alias for resize()
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 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.