From 03f0092cefadc86757e89bc1da9f2be678591d3d Mon Sep 17 00:00:00 2001 From: Henry <Henry> Date: Mon, 1 Nov 2010 11:07:08 +0000 Subject: [PATCH] fvcSmooth: The alpha field discriminant parameters are now optional arguments --- .../finiteVolume/fvc/fvcSmooth/fvcSmooth.C | 27 ++++++++++++------- .../finiteVolume/fvc/fvcSmooth/fvcSmooth.H | 23 ++++++++++------ 2 files changed, 32 insertions(+), 18 deletions(-) diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C index 47cfb5582ac..4ef505b3f90 100644 --- a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C +++ b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C @@ -123,7 +123,10 @@ void Foam::fvc::spread ( volScalarField& field, const volScalarField& alpha, - const label nLayers + const label nLayers, + const scalar alphaDiff, + const scalar alphaMax, + const scalar alphaMin ) { const fvMesh& mesh = field.mesh(); @@ -153,11 +156,11 @@ void Foam::fvc::spread if ( - (alpha[own] > 0.01 && alpha[own] < 0.99) - || (alpha[nbr] > 0.01 && alpha[nbr] < 0.99) + (alpha[own] > alphaMin && alpha[own] < alphaMax) + || (alpha[nbr] > alphaMin && alpha[nbr] < alphaMax) ) { - if (mag(alpha[own] - alpha[nbr]) > 0.2) + if (mag(alpha[own] - alpha[nbr]) > alphaDiff) { changedFaces.append(facei); changedFacesInfo.append @@ -185,11 +188,14 @@ void Foam::fvc::spread if ( - (alpha[own] > 0.01 && alpha[own] < 0.99) - || (alphapn[patchFacei] > 0.01 && alphapn[patchFacei] < 0.99) + (alpha[own] > alphaMin && alpha[own] < alphaMax) + || ( + alphapn[patchFacei] > alphaMin + && alphapn[patchFacei] < alphaMax + ) ) { - if (mag(alpha[own] - alphapn[patchFacei]) > 0.2) + if (mag(alpha[own] - alphapn[patchFacei]) > alphaDiff) { changedFaces.append(facei); changedFacesInfo.append(smoothData(field[own])); @@ -227,7 +233,8 @@ void Foam::fvc::sweep ( volScalarField& field, const volScalarField& alpha, - const label nLayers + const label nLayers, + const scalar alphaDiff ) { const fvMesh& mesh = field.mesh(); @@ -250,7 +257,7 @@ void Foam::fvc::sweep const label own = owner[facei]; const label nbr = neighbour[facei]; - if (mag(alpha[own] - alpha[nbr]) > 0.2) + if (mag(alpha[own] - alpha[nbr]) > alphaDiff) { changedFaces.append(facei); changedFacesInfo.append @@ -275,7 +282,7 @@ void Foam::fvc::sweep scalarField alphapn = alpha.boundaryField()[patchi].patchNeighbourField(); - if (mag(alpha[own] - alphapn[patchFacei]) > 0.2) + if (mag(alpha[own] - alphapn[patchFacei]) > alphaDiff) { changedFaces.append(facei); changedFacesInfo.append diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.H b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.H index 2ce5b1b5d32..ef65494c459 100644 --- a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.H +++ b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.H @@ -28,15 +28,18 @@ Description Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistribute the first field argument. - smooth: smooths the field by ensuring the values in neighbouring - cells are at least coeff* the cell value. + smooth: smooths the field by ensuring the values in neighbouring cells are + at least coeff* the cell value. - spread: redistributes the field by spreading the maximum value within - the region defined by the value and gradient of alpha. + spread: redistributes the field by spreading the maximum value within the + region defined by the value (being between alphaMax and alphaMin) + and gradient of alpha (where the difference between the values in + neighbouring cells is larger than alphaDiff). sweep: redistributes the field by sweeping the maximum value where the - gradient of alpha is large away from that starting point of the - sweep. + gradient of alpha is large (where the difference between the values + in neighbouring cells is larger than alphaDiff) away from that + starting point of the sweep. SourceFiles fvcSmooth.C @@ -64,14 +67,18 @@ namespace fvc ( volScalarField& field, const volScalarField& alpha, - const label nLayers + const label nLayers, + const scalar alphaDiff = 0.2, + const scalar alphaMax = 0.99, + const scalar alphaMin = 0.01 ); void sweep ( volScalarField& field, const volScalarField& alpha, - const label nLayers + const label nLayers, + const scalar alphaDiff = 0.2 ); } } -- GitLab