Go to the documentation of this file.
35 template<
class Polynomial>
40 const scalar linearLimitFactor,
41 const scalar centralWeight
51 mesh, stencil,
true, linearLimitFactor, centralWeight
53 coeffs_(
mesh.nFaces())
58 <<
"Contructing CentredFitSnGradData<Polynomial>" <<
endl;
65 Info<<
" Finished constructing polynomialFit data" <<
endl;
72 template<
class Polynomial>
78 const scalar deltaCoeff,
85 this->findFaceDirs(idir, jdir, kdir, facei);
89 wts[0] = this->centralWeight();
90 wts[1] = this->centralWeight();
121 Polynomial::addCoeffs(
B[ip], d, wts[ip], this->dim());
125 for (
label i = 0; i <
B.m(); i++)
132 label stencilSize =
C.size();
135 bool goodFit =
false;
136 for (
int iIt = 0; iIt < 8 && !goodFit; iIt++)
141 for (
label i=0; i<stencilSize; i++)
143 coeffsi[i] = wts[1]*wts[i]*invB(1, i)/scale;
148 mag(wts[0]*wts[0]*invB(0, 0) - wLin)
149 < this->linearLimitFactor()*wLin)
150 && (
mag(wts[0]*wts[1]*invB(0, 1) - (1 - wLin)
151 ) < this->linearLimitFactor()*(1 - wLin))
152 && coeffsi[0] < 0 && coeffsi[1] > 0
153 &&
mag(coeffsi[0] + deltaCoeff) < 0.5*deltaCoeff
154 &&
mag(coeffsi[1] - deltaCoeff) < 0.5*deltaCoeff;
162 <<
"Cannot fit face " << facei <<
" iteration " << iIt
163 <<
" with sum of weights " <<
sum(coeffsi) <<
nl
164 <<
" Weights " << coeffsi <<
nl
165 <<
" Linear weights " << wLin <<
" " << 1 - wLin <<
nl
166 <<
" deltaCoeff " << deltaCoeff <<
nl
167 <<
" sing vals " << svd.
S() <<
nl
168 <<
"Components of goodFit:\n"
169 <<
" wts[0]*wts[0]*invB(0, 0) = "
170 << wts[0]*wts[0]*invB(0, 0) <<
nl
171 <<
" wts[0]*wts[1]*invB(0, 1) = "
172 << wts[0]*wts[1]*invB(0, 1)
173 <<
" dim = " << this->dim() <<
endl;
178 for (
label j = 0; j <
B.n(); j++)
184 for (
label i = 0; i <
B.m(); i++)
195 coeffsi[0] += deltaCoeff;
196 coeffsi[1] -= deltaCoeff;
201 <<
"Could not fit face " << facei
202 <<
" Coefficients = " << coeffsi
203 <<
", reverting to uncorrected." <<
endl;
210 template<
class Polynomial>
219 this->stencil().collectData(
mesh.C(), stencilPoints);
226 for (
label facei = 0; facei <
mesh.nInternalFaces(); facei++)
231 stencilPoints[facei],
239 const surfaceScalarField::Boundary& bdC = dC.
boundaryField();
255 stencilPoints[facei],
int debug
Static debugging option.
const Cmpt & x() const
Access to the vector x component.
#define InfoInFunction
Report an information message using Foam::Info.
Data for the upwinded and centred polynomial fit interpolation schemes. The linearCorrection_ determi...
static constexpr const zero Zero
Global zero.
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)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
scalarRectangularMatrix VSinvUt() const
Return the matrix product V S^(-1) U^T (the pseudo inverse)
messageStream Info
Information stream (uses stdout - output is on the master only)
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):
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.