Commit 1c101ec1 authored by Mark Olesen's avatar Mark Olesen
Browse files

bugfix fixedBlended scheme, cosmetics changes

- actually blend correction factors as well in fixedBlended scheme

- consistency change: move 'linearUpwindV' into same directory as
  'linearUpwind', move correction code from .H into .C code.

- minor docu format/spelling changes

A more general question:
~~~~~~~~~~~~~~~~~~~~~~~~
- would it make more sense to use template specializations for the
  '*V' differencing schemes?
parent d2a78f7e
......@@ -206,7 +206,7 @@ $(limitedSchemes)/limitedSurfaceInterpolationScheme/limitedSurfaceInterpolationS
$(limitedSchemes)/upwind/upwind.C
$(limitedSchemes)/blended/blended.C
$(limitedSchemes)/linearUpwind/linearUpwind.C
$(limitedSchemes)/linearUpwindV/linearUpwindV.C
$(limitedSchemes)/linearUpwind/linearUpwindV.C
$(limitedSchemes)/Gamma/Gamma.C
$(limitedSchemes)/SFCD/SFCD.C
$(limitedSchemes)/Minmod/Minmod.C
......
......@@ -21,7 +21,7 @@ License
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "cellLimitedGrad.H"
......@@ -377,7 +377,7 @@ tmp<volTensorField> cellLimitedGrad<vector>::grad
);
}
}
if (fv::debug)
{
Info<< "gradient limiter for: " << vsf.name()
......
......@@ -21,7 +21,7 @@ License
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "cellMDLimitedGrad.H"
......
......@@ -41,7 +41,7 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
Class LimitedLimiter Declaration
Class LimitedLimiter Declaration
\*---------------------------------------------------------------------------*/
template<class LimitedScheme>
......
......@@ -41,7 +41,7 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
Class Limited01Limiter Declaration
Class Limited01Limiter Declaration
\*---------------------------------------------------------------------------*/
template<class LimitedScheme>
......
......@@ -79,18 +79,17 @@ public:
{
scalar phiCD = cdWeight*phiP + (1 - cdWeight)*phiN;
scalar phif;
scalar phiU;
scalar phiU, phif;
if (faceFlux > 0)
{
phif = 0.5*(phiCD + phiP + (1 - cdWeight)*(d & gradcP));
phiU = phiP;
phif = 0.5*(phiCD + phiP + (1 - cdWeight)*(d & gradcP));
}
else
{
phif = 0.5*(phiCD + phiN - cdWeight*(d & gradcN));
phiU = phiN;
phif = 0.5*(phiCD + phiN - cdWeight*(d & gradcN));
}
// Calculate the effective limiter for the QUICK interpolation
......
......@@ -81,8 +81,7 @@ public:
scalar phiCD = gradfV & (cdWeight*phiP + (1 - cdWeight)*phiN);
scalar phif;
scalar phiU;
scalar phiU, phif;
if (faceFlux > 0)
{
......@@ -95,7 +94,7 @@ public:
phif = 0.5*(phiCD + phiU - cdWeight*(gradfV & (d & gradcN)));
}
// Calculate the effective limiter for the linearUpwind interpolation
// Calculate the effective limiter for the QUICK interpolation
scalar QLimiter = (phif - phiU)/stabilise(phiCD - phiU, SMALL);
// Limit the limiter between upwind and downwind
......
......@@ -26,7 +26,7 @@ Class
Foam::filteredLinearLimiter
Description
Class to generate weighting factors for the filtered-linear
Class to generate weighting factors for the filteredLinear
differencing scheme.
The aim is to remove high-frequency modes with "staggering"
......
......@@ -26,7 +26,7 @@ Class
Foam::filteredLinear2Limiter
Description
Class to generate weighting factors for the filtered-linear-2
Class to generate weighting factors for the filteredLinear2
differencing scheme.
The aim is to remove high-frequency modes with "staggering"
......
......@@ -26,11 +26,13 @@ Class
Foam::filteredLinear2VLimiter
Description
Class to generate weighting factors for the filteredLinear2V differencing
scheme. The aim is to remove high-frequency modes with "staggering"
characteristics from vector fields by comparing the face gradient in the
direction of maximum gradient with both neighbouring cell gradients and
introduce small amounts of upwind in order to damp these modes.
Class to generate weighting factors for the filteredLinear2V
differencing scheme.
The aim is to remove high-frequency modes with "staggering"
characteristics from vector fields by comparing the face gradient in
the direction of maximum gradient with both neighbouring cell gradients
and introduce small amounts of upwind in order to damp these modes.
Used in conjunction with the template class LimitedScheme.
......
......@@ -26,7 +26,7 @@ Class
Foam::filteredLinear3Limiter
Description
Class to generate weighting factors for the filtered-linear-3
Class to generate weighting factors for the filteredLinear
differencing scheme.
The aim is to remove high-frequency modes with "staggering"
......@@ -79,7 +79,7 @@ public:
<< "coefficient = " << k_
<< " should be >= 0 and <= 1"
<< exit(FatalIOError);
}
}
}
scalar limiter
......
......@@ -26,11 +26,13 @@ Class
Foam::filteredLinear3VLimiter
Description
Class to generate weighting factors for the filteredLinear3V differencing
scheme. The aim is to remove high-frequency modes with "staggering"
characteristics from vector fields by comparing the face gradient in the
direction of maximum gradient with both neighbouring cell gradients and
introduce small amounts of upwind in order to damp these modes.
Class to generate weighting factors for the filteredLinear3V
differencing scheme.
The aim is to remove high-frequency modes with "staggering"
characteristics from vector fields by comparing the face gradient in
the direction of maximum gradient with both neighbouring cell gradients
and introduce small amounts of upwind in order to damp these modes.
Used in conjunction with the template class LimitedScheme.
......@@ -77,7 +79,7 @@ public:
<< "coefficient = " << k_
<< " should be >= 0 and <= 1"
<< exit(FatalIOError);
}
}
}
scalar limiter
......
......@@ -90,7 +90,7 @@ public:
)
{}
//- Construct from Istream.
//- Construct from Istream.
// The name of the flux field is read from the Istream and looked-up
// from the mesh objectRegistry
linearUpwind
......@@ -142,8 +142,9 @@ public:
virtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
correction
(
const GeometricField<Type, fvPatchField, volMesh>& vf
const GeometricField<Type, fvPatchField, volMesh>&
) const;
};
......
......@@ -31,6 +31,103 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh> >
Foam::linearUpwindV<Type>::correction
(
const GeometricField<Type, fvPatchField, volMesh>& vf
) const
{
const fvMesh& mesh = this->mesh();
tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsfCorr
(
new GeometricField<Type, fvsPatchField, surfaceMesh>
(
IOobject
(
vf.name(),
mesh.time().timeName(),
mesh
),
mesh,
dimensioned<Type>
(
vf.name(),
vf.dimensions(),
pTraits<Type>::zero
)
)
);
GeometricField<Type, fvsPatchField, surfaceMesh>& sfCorr = tsfCorr();
const surfaceScalarField& faceFlux = this->faceFlux_;
const surfaceScalarField& w = mesh.weights();
const labelList& own = mesh.owner();
const labelList& nei = mesh.neighbour();
const vectorField& C = mesh.C();
const vectorField& Cf = mesh.Cf();
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
gradVf = gradScheme_().grad(vf);
forAll(faceFlux, facei)
{
vector maxCorr;
if (faceFlux[facei] > 0.0)
{
maxCorr =
(1.0 - w[facei])
*(vf[nei[facei]] - vf[own[facei]]);
sfCorr[facei] =
(Cf[facei] - C[own[facei]]) & gradVf[own[facei]];
}
else
{
maxCorr =
w[facei]*(vf[own[facei]] - vf[nei[facei]]);
sfCorr[facei] =
(Cf[facei] - C[nei[facei]]) & gradVf[nei[facei]];
}
scalar sfCorrs = magSqr(sfCorr[facei]);
scalar maxCorrs = sfCorr[facei] & maxCorr;
if (sfCorrs > 0)
{
if (maxCorrs < 0)
{
sfCorr[facei] = vector::zero;
}
else if (sfCorrs > maxCorrs)
{
sfCorr[facei] *= maxCorrs/(sfCorrs + VSMALL);
}
}
else if (sfCorrs < 0)
{
if (maxCorrs > 0)
{
sfCorr[facei] = vector::zero;
}
else if (sfCorrs < maxCorrs)
{
sfCorr[facei] *= maxCorrs/(sfCorrs - VSMALL);
}
}
}
return tsfCorr;
}
namespace Foam
{
makelimitedSurfaceInterpolationTypeScheme(linearUpwindV, vector)
......
......@@ -90,7 +90,7 @@ public:
)
{}
//- Construct from Istream.
//- Construct from Istream.
// The name of the flux field is read from the Istream and looked-up
// from the mesh objectRegistry
linearUpwindV
......@@ -142,97 +142,9 @@ public:
virtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
correction
(
const GeometricField<Type, fvPatchField, volMesh>& vf
) const
{
const fvMesh& mesh = this->mesh();
const GeometricField<Type, fvPatchField, volMesh>&
) const;
tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsfCorr
(
new GeometricField<Type, fvsPatchField, surfaceMesh>
(
IOobject
(
vf.name(),
mesh.time().timeName(),
mesh
),
mesh,
dimensioned<Type>
(
vf.name(),
vf.dimensions(),
pTraits<Type>::zero
)
)
);
GeometricField<Type, fvsPatchField, surfaceMesh>& sfCorr = tsfCorr();
const surfaceScalarField& faceFlux = this->faceFlux_;
const surfaceScalarField& w = mesh.weights();
const labelList& own = mesh.owner();
const labelList& nei = mesh.neighbour();
const vectorField& C = mesh.C();
const vectorField& Cf = mesh.Cf();
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
gradVf = gradScheme_().grad(vf);
forAll(faceFlux, facei)
{
vector maxCorr;
if (faceFlux[facei] > 0.0)
{
maxCorr =
(1.0 - w[facei])
*(vf[nei[facei]] - vf[own[facei]]);
sfCorr[facei] =
(Cf[facei] - C[own[facei]]) & gradVf[own[facei]];
}
else
{
maxCorr =
w[facei]*(vf[own[facei]] - vf[nei[facei]]);
sfCorr[facei] =
(Cf[facei] - C[nei[facei]]) & gradVf[nei[facei]];
}
scalar sfCorrs = magSqr(sfCorr[facei]);
scalar maxCorrs = sfCorr[facei] & maxCorr;
if (sfCorrs > 0)
{
if (maxCorrs < 0)
{
sfCorr[facei] = vector::zero;
}
else if (sfCorrs > maxCorrs)
{
sfCorr[facei] *= maxCorrs/(sfCorrs + VSMALL);
}
}
else if (sfCorrs < 0)
{
if (maxCorrs > 0)
{
sfCorr[facei] = vector::zero;
}
else if (sfCorrs < maxCorrs)
{
sfCorr[facei] *= maxCorrs/(sfCorrs - VSMALL);
}
}
}
return tsfCorr;
}
};
......
......@@ -58,11 +58,11 @@ class UpwindFitData
// Private data
//- For each face of the mesh store the coefficients to multiply the
// stencil cell values by if the flow is from the owner
// stencil cell values by if the flow is from the owner
List<scalarList> owncoeffs_;
//- For each face of the mesh store the coefficients to multiply the
// stencil cell values by if the flow is from the neighbour
// stencil cell values by if the flow is from the neighbour
List<scalarList> neicoeffs_;
......
......@@ -26,7 +26,7 @@ Class
Foam::UpwindFitScheme
Description
Upwind biased fit surface interpolation scheme which applies an explicit
Upwind biased fit surface interpolation scheme that applies an explicit
correction to linear.
\*---------------------------------------------------------------------------*/
......
......@@ -27,7 +27,8 @@ Class
Description
BiLinear polynomial for interpolation fitting.
Can be used with the CentredFit scheme to crate a biLinear surface
Can be used with the CentredFit scheme to create a biLinear surface
interpolation scheme
\*---------------------------------------------------------------------------*/
......
......@@ -28,7 +28,7 @@ Class
Description
Cubic polynomial for upwind biased interpolation fitting.
Can be used with the UpwindFit scheme to crate a cubic surface
Can be used with the UpwindFit scheme to create a cubic surface
interpolation scheme
\*---------------------------------------------------------------------------*/
......
......@@ -125,6 +125,7 @@ public:
}
}
//- Construct from mesh, faceFlux and Istream
fixedBlended
(
......@@ -165,7 +166,8 @@ public:
// Member Functions
//- Return the interpolation weighting factors
tmp<surfaceScalarField> weights
tmp<surfaceScalarField>
weights
(
const GeometricField<Type, fvPatchField, volMesh>& vf
) const
......@@ -175,15 +177,74 @@ public:
+ (scalar(1.0) - blendingFactor_)*tScheme2_().weights(vf);
}
//- Return the face-interpolate of the given cell field
// with explicit correction
tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
interpolate(const GeometricField<Type, fvPatchField, volMesh>& vf) const
interpolate
(
const GeometricField<Type, fvPatchField, volMesh>& vf
) const
{
return
blendingFactor_*tScheme1_().interpolate(vf)
+ (scalar(1.0) - blendingFactor_)*tScheme2_().interpolate(vf);
}
//- Return true if this scheme uses an explicit correction
virtual bool corrected() const
{
return tScheme1_().corrected() || tScheme2_().corrected();
}
//- Return the explicit correction to the face-interpolate
// for the given field
virtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
correction
(
const GeometricField<Type, fvPatchField, volMesh>& vf
) const
{
if (tScheme1_().corrected())
{
if (tScheme2_().corrected())
{
return
(
blendingFactor_
* tScheme1_().correction(vf)
+ (scalar(1.0) - blendingFactor_)
* tScheme2_().correction(vf)
);
}
else
{
return
(
blendingFactor_
* tScheme1_().correction(vf)
);
}
}
else if (tScheme2_().corrected())
{
return
(
(scalar(1.0) - blendingFactor_)
* tScheme2_().correction(vf)
);
}
else
{
return tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
(
NULL
);
}
}
};
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment