[MichiB]

Hi,

I have had a look at the implicit discretization of anisotropic diffusion (fvm::laplacian).
It seems to me that the implementation is wrong. At least I haven't yet managed to figure out why the coefficients are as they are.
It is fairly easy to follow in Ferzinger & Peric how the non-orthogonal correction approach is done for scalar "gamma" fields.
For tensor fields it's a bit more tricky.
What is being passed to the uncorrected (implicit) laplacian function is the scalar field sfGammaSn.
For me it would make sense to pass a vector field sfGamma and then split the derivative nabla(vf) into a derivative along d and a derivative along
sfGamma_{unity}.

I think the way fvm::laplacian is implemented now shouldn't even work for scalar "gammas" right?
Otherwise there wouldn't exist a template version for fvm::laplacians for scalars.

I would appriciate if anybody can exchange thoughts on this topic / approves my thoughts / teaches me wrong / explains to me the implementation.

Cheers,
Michi

template<class Type, class GType>
tmp<fvMatrix<Type> >
gaussLaplacianScheme<Type, GType>::fvmLaplacianUncorrected
(
const surfaceScalarField& gammaMagSf, actually passes sfGammaSn
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
tmp<surfaceScalarField> tdeltaCoeffs =
this->tsnGradScheme_().deltaCoeffs(vf);
const surfaceScalarField& deltaCoeffs = tdeltaCoeffs();

tmp<fvMatrix<Type> > tfvm
(
new fvMatrix<Type>
(
vf,
deltaCoeffs.dimensions()*gammaMagSf.dimensions()*v f.dimensions()
)
);
fvMatrix<Type>& fvm = tfvm();

fvm.upper() = deltaCoeffs.internalField()*gammaMagSf.internalFie ld();
fvm.negSumDiag();

[MichiB]

I figured out that the coeffs indeed correspond to what I've derived

which is (Sf dot Gamma dot Sf)/(Sf dot d) ... for the uncorrected laplacian.

I hope I didn't confuse anyone that read this message! If so, sorry!