Added flux handling to the correction functions.

This commit is contained in:
henry
2009-03-30 13:20:56 +01:00
parent 002c85dd67
commit fc33d49ae0
2 changed files with 44 additions and 31 deletions

View File

@ -57,22 +57,7 @@ Foam::tmp<Foam::Field<Type> > Foam::lduMatrix::H(const Field<Type>& psi) const
for (register label face=0; face<nFaces; face++)
{
#ifdef ICC_IA64_PREFETCH
__builtin_prefetch (&uPtr[face+32],0,0);
__builtin_prefetch (&lPtr[face+32],0,0);
__builtin_prefetch (&lowerPtr[face+32],0,1);
__builtin_prefetch (&psiPtr[lPtr[face+32]],0,1);
__builtin_prefetch (&HpsiPtr[uPtr[face+32]],0,1);
#endif
HpsiPtr[uPtr[face]] -= lowerPtr[face]*psiPtr[lPtr[face]];
#ifdef ICC_IA64_PREFETCH
__builtin_prefetch (&upperPtr[face+32],0,1);
__builtin_prefetch (&psiPtr[uPtr[face+32]],0,1);
__builtin_prefetch (&HpsiPtr[lPtr[face+32]],0,1);
#endif
HpsiPtr[lPtr[face]] -= upperPtr[face]*psiPtr[uPtr[face]];
}
}
@ -94,22 +79,33 @@ template<class Type>
Foam::tmp<Foam::Field<Type> >
Foam::lduMatrix::faceH(const Field<Type>& psi) const
{
const scalarField& Lower = const_cast<const lduMatrix&>(*this).lower();
const scalarField& Upper = const_cast<const lduMatrix&>(*this).upper();
// Take refereces to addressing
const unallocLabelList& l = lduAddr().lowerAddr();
const unallocLabelList& u = lduAddr().upperAddr();
tmp<Field<Type> > tfaceHpsi(new Field<Type> (Lower.size()));
Field<Type> & faceHpsi = tfaceHpsi();
for (register label face=0; face<l.size(); face++)
if (lowerPtr_ || upperPtr_)
{
faceHpsi[face] = Upper[face]*psi[u[face]] - Lower[face]*psi[l[face]];
}
const scalarField& Lower = const_cast<const lduMatrix&>(*this).lower();
const scalarField& Upper = const_cast<const lduMatrix&>(*this).upper();
return tfaceHpsi;
const unallocLabelList& l = lduAddr().lowerAddr();
const unallocLabelList& u = lduAddr().upperAddr();
tmp<Field<Type> > tfaceHpsi(new Field<Type> (Lower.size()));
Field<Type> & faceHpsi = tfaceHpsi();
for (register label face=0; face<l.size(); face++)
{
faceHpsi[face] = Upper[face]*psi[u[face]] - Lower[face]*psi[l[face]];
}
return tfaceHpsi;
}
else
{
FatalErrorIn("lduMatrix::faceH(const Field<Type>& psi) const")
<< "Cannot calculate faceH"
" the matrix does not have any off-diagonal coefficients."
<< exit(FatalError);
return tmp<Field<Type> >(NULL);
}
}

View File

@ -1303,7 +1303,14 @@ Foam::tmp<Foam::fvMatrix<Type> > Foam::correction
const fvMatrix<Type>& A
)
{
return A - (A & A.psi());
tmp<Foam::fvMatrix<Type> > tAcorr = A - (A & A.psi());
if ((A.hasUpper() || A.hasLower()) && A.mesh().fluxRequired(A.psi().name()))
{
tAcorr().faceFluxCorrectionPtr() = (-A.flux()).ptr();
}
return tAcorr;
}
@ -1313,7 +1320,17 @@ Foam::tmp<Foam::fvMatrix<Type> > Foam::correction
const tmp<fvMatrix<Type> >& tA
)
{
return tA - (tA() & tA().psi());
tmp<Foam::fvMatrix<Type> > tAcorr = tA - (tA() & tA().psi());
// Note the matrix coefficients are still that of matrix A
const fvMatrix<Type>& A = tAcorr();
if ((A.hasUpper() || A.hasLower()) && A.mesh().fluxRequired(A.psi().name()))
{
tAcorr().faceFluxCorrectionPtr() = (-A.flux()).ptr();
}
return tAcorr;
}