diff --git a/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C b/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C index 5f25016a84eeed80d3545cc386be4b50842ddbbe..a4282534217c6742fd2739cea7645961cf478f2d 100644 --- a/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C +++ b/applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -165,7 +165,6 @@ int main(int argc, char *argv[]) mesh.movePoints(motionPtr->newPoints()); phi = aphiv_pos*rho_pos + aphiv_neg*rho_neg; - Info<< phi.boundaryField()[0] << endl; surfaceVectorField phiUp ( diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearUpwind/linearUpwindV.C b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearUpwind/linearUpwindV.C index 9e5ad7f3a02e89bb57288775d34f7ddc3ee46f97..6368d2accd863cc9d6965b94b0a826e6ef2bbbaa 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearUpwind/linearUpwindV.C +++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearUpwind/linearUpwindV.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -70,8 +70,8 @@ Foam::linearUpwindV<Type>::correction const labelList& own = mesh.owner(); const labelList& nei = mesh.neighbour(); - const vectorField& C = mesh.C(); - const vectorField& Cf = mesh.Cf(); + const volVectorField& C = mesh.C(); + const surfaceVectorField& Cf = mesh.Cf(); tmp < @@ -97,8 +97,7 @@ Foam::linearUpwindV<Type>::correction if (faceFlux[facei] > 0.0) { maxCorr = - (1.0 - w[facei]) - *(vf[nei[facei]] - vf[own[facei]]); + (1.0 - w[facei])*(vf[nei[facei]] - vf[own[facei]]); sfCorr[facei] = (Cf[facei] - C[own[facei]]) & gradVf[own[facei]]; @@ -139,6 +138,86 @@ Foam::linearUpwindV<Type>::correction } } + + typename GeometricField<Type, fvsPatchField, surfaceMesh>:: + GeometricBoundaryField& bSfCorr = sfCorr.boundaryField(); + + forAll(bSfCorr, patchi) + { + fvsPatchField<Type>& pSfCorr = bSfCorr[patchi]; + + if (pSfCorr.coupled()) + { + const labelUList& pOwner = + mesh.boundary()[patchi].faceCells(); + + const vectorField& pCf = Cf.boundaryField()[patchi]; + const scalarField& pW = w.boundaryField()[patchi]; + + const scalarField& pFaceFlux = faceFlux.boundaryField()[patchi]; + + const Field<typename outerProduct<vector, Type>::type> pGradVfNei + ( + gradVf.boundaryField()[patchi].patchNeighbourField() + ); + + const Field<Type> pVfNei + ( + vf.boundaryField()[patchi].patchNeighbourField() + ); + + // Build the d-vectors + vectorField pd(Cf.boundaryField()[patchi].patch().delta()); + + forAll(pOwner, facei) + { + label own = pOwner[facei]; + + vector maxCorr; + + if (pFaceFlux[facei] > 0) + { + pSfCorr[facei] = (pCf[facei] - C[own]) & gradVf[own]; + + maxCorr = (1.0 - pW[facei])*(pVfNei[facei] - vf[own]); + } + else + { + pSfCorr[facei] = + (pCf[facei] - pd[facei] - C[own]) & pGradVfNei[facei]; + + maxCorr = pW[facei]*(vf[own] - pVfNei[facei]); + } + + scalar pSfCorrs = magSqr(pSfCorr[facei]); + scalar maxCorrs = pSfCorr[facei] & maxCorr; + + if (pSfCorrs > 0) + { + if (maxCorrs < 0) + { + pSfCorr[facei] = vector::zero; + } + else if (pSfCorrs > maxCorrs) + { + pSfCorr[facei] *= maxCorrs/(pSfCorrs + VSMALL); + } + } + else if (pSfCorrs < 0) + { + if (maxCorrs > 0) + { + pSfCorr[facei] = vector::zero; + } + else if (pSfCorrs < maxCorrs) + { + pSfCorr[facei] *= maxCorrs/(pSfCorrs - VSMALL); + } + } + } + } + } + return tsfCorr; }