Bug fix: harmonic interpolation @ coupled BC / optimisation

This commit is contained in:
Henrik Rusche 2011-03-18 12:25:39 +01:00
parent 4a271dbf06
commit 7a5c10e48d
2 changed files with 52 additions and 20 deletions

View file

@ -180,25 +180,27 @@ void regionCoupleFvPatchField<Type>::evaluate
Field<Type> fOwn = this->patchInternalField();
Field<Type> fNei = this->patchNeighbourField();
scalarField mOwn = mag(fOwn);
scalarField mNei = mag(fNei);
scalarField mean = 2*(mOwn*mNei)/(mOwn + mNei + SMALL);
scalarField magFOwn = mag(fOwn);
scalarField magFNei = mag(fNei);
scalarField weights(fOwn.size(), 0.5);
scalar den;
// Calculate internal weights using field magnitude
scalarField weights(fOwn.size());
forAll (weights, faceI)
{
den = mOwn[faceI] - mNei[faceI];
scalar mOwn = magFOwn[faceI];
scalar mNei = magFNei[faceI];
scalar den = mOwn - mNei;
if (mag(den) > SMALL)
{
weights[faceI] = (mean[faceI] - mNei[faceI])/den;
scalar mean = 2.0*mOwn*mNei/(mOwn + mNei);
weights[faceI] = (mean - mNei)/den;
}
else
{
// Use 0.5 weights
weights[faceI] = 0.5;
}
}

View file

@ -131,33 +131,63 @@ public:
scalarField magPhi = mag(phi);
// Initialise weights to 0.5 for uniform field (den = 0)
scalarField& wIn = w.internalField();
wIn = 0.5;
// Calculate internal weights using field magnitude
scalar mOwn, mNei, den, mean;
forAll (owner, faceI)
{
mOwn = magPhi[owner[faceI]];
mNei = magPhi[neighbour[faceI]];
scalar mOwn = magPhi[owner[faceI]];
scalar mNei = magPhi[neighbour[faceI]];
mean = 2*(mOwn*mNei)/(mOwn + mNei + SMALL);
den = mOwn - mNei;
scalar den = mOwn - mNei;
if (mag(den) > SMALL)
{
scalar mean = 2.0*mOwn*mNei/(mOwn + mNei);
wIn[faceI] = (mean - mNei)/den;
}
else
{
// Use 0.5 weights
wIn[faceI] = 0.5;
}
}
// Boundary weights are 1
w.boundaryField() = 1;
forAll (phi.boundaryField(), pi)
{
fvsPatchScalarField& wp = w.boundaryField()[pi];
const fvPatchField<Type>& patchPhi = phi.boundaryField()[pi];
if (patchPhi.coupled())
{
scalarField magPhiOwn = mag(patchPhi.patchInternalField());
scalarField magPhiNei = mag(patchPhi.patchNeighbourField());
// Calculate internal weights using field magnitude
forAll (patchPhi, faceI)
{
scalar mOwn = magPhiOwn[faceI];
scalar mNei = magPhiNei[faceI];
scalar den = mOwn - mNei;
if (mag(den) > SMALL)
{
scalar mean = 2.0*mOwn*mNei/(mOwn + mNei);
wp[faceI] = (mean - mNei)/den;
}
else
{
wp[faceI] = 0.5;
}
}
}
else
{
// Boundary weights for uncoupled patches are 1
wp = 1;
}
}
return tw;
}