Commit 592692a0 authored by Mark Olesen's avatar Mark Olesen
Browse files

ENH: added surfaceFieldValue uniformity operation

parent dff69e39
......@@ -3,7 +3,7 @@
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation | Copyright (C) 2017 OpenCFD Ltd.
\\/ M anipulation | Copyright (C) 2017-2018 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
......@@ -83,18 +83,21 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::operationTypeNames_
{ operationType::opCoV, "CoV" },
{ operationType::opAreaNormalAverage, "areaNormalAverage" },
{ operationType::opAreaNormalIntegrate, "areaNormalIntegrate" },
{ operationType::opUniformity, "uniformity" },
// Using weighting
{ operationType::opWeightedSum, "weightedSum" },
{ operationType::opWeightedAverage, "weightedAverage" },
{ operationType::opWeightedAreaAverage, "weightedAreaAverage" },
{ operationType::opWeightedAreaIntegrate, "weightedAreaIntegrate" },
{ operationType::opWeightedUniformity, "weightedUniformity" },
// Using absolute weighting
{ operationType::opAbsWeightedSum, "absWeightedSum" },
{ operationType::opAbsWeightedAverage, "absWeightedAverage" },
{ operationType::opAbsWeightedAreaAverage, "absWeightedAreaAverage" },
{ operationType::opAbsWeightedAreaIntegrate, "absWeightedAreaIntegrate" },
{ operationType::opAbsWeightedUniformity, "absWeightedUniformity" },
};
const Foam::Enum
......@@ -117,10 +120,8 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::obr() const
{
return mesh_.lookupObject<objectRegistry>(regionName_);
}
else
{
return mesh_;
}
return mesh_;
}
......@@ -599,7 +600,7 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::initialise
}
}
// Backwards compatibility for v1612+ and older
// Backwards compatibility for v1612 and older
List<word> orientedFields;
if (dict.readIfPresent("orientedFields", orientedFields))
{
......@@ -696,6 +697,45 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::processValues
return gSum(pos0(nv)*mag(values) - neg(nv)*mag(values));
}
case opUniformity:
case opWeightedUniformity:
case opAbsWeightedUniformity:
{
const scalar areaTotal = gSum(mag(Sf));
tmp<scalarField> areaVal(values * mag(Sf));
scalar mean, numer;
if (canWeight(weightField))
{
// Weighted quantity = (Weight * phi * dA)
tmp<scalarField> weight(weightingFactor(weightField));
// Mean weighted value (area-averaged)
mean = gSum(weight()*areaVal()) / areaTotal;
// Abs. deviation from weighted mean value
numer = gSum(mag(weight*areaVal - (mean * mag(Sf))));
}
else
{
// Unweighted quantity = (1 * phi * dA)
// Mean value (area-averaged)
mean = gSum(areaVal()) / areaTotal;
// Abs. deviation from mean value
numer = gSum(mag(areaVal - (mean * mag(Sf))));
}
// Uniformity index
const scalar ui = 1 - numer/(2*mag(mean*areaTotal) + ROOTVSMALL);
return min(max(ui, 0), 1);
}
default:
{
// Fall through to other operations
......@@ -742,6 +782,45 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::processValues
const scalar val = gSum(values & Sf);
return vector(val, 0, 0);
}
case opUniformity:
case opWeightedUniformity:
case opAbsWeightedUniformity:
{
const scalar areaTotal = gSum(mag(Sf));
tmp<scalarField> areaVal(values & Sf);
scalar mean, numer;
if (canWeight(weightField))
{
// Weighted quantity = (Weight * phi . dA)
tmp<scalarField> weight(weightingFactor(weightField));
// Mean weighted value (area-averaged)
mean = gSum(weight()*areaVal()) / areaTotal;
// Abs. deviation from weighted mean value
numer = gSum(mag(weight*areaVal - (mean * mag(Sf))));
}
else
{
// Unweighted quantity = (1 * phi . dA)
// Mean value (area-averaged)
mean = gSum(areaVal()) / areaTotal;
// Abs. deviation from mean value
numer = gSum(mag(areaVal - (mean * mag(Sf))));
}
// Uniformity index
const scalar ui = 1 - numer/(2*mag(mean*areaTotal) + ROOTVSMALL);
return vector(min(max(ui, 0), 1), 0, 0);
}
default:
{
// Fall through to other operations
......@@ -762,11 +841,9 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
{
return mag(weightField);
}
else
{
// pass through
return weightField;
}
// pass through
return weightField;
}
......@@ -788,10 +865,8 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
{
return mag(weightField * mag(Sf));
}
else
{
return (weightField * mag(Sf));
}
return (weightField * mag(Sf));
}
......@@ -813,10 +888,8 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::weightingFactor
{
return mag(weightField & Sf);
}
else
{
return (weightField & Sf);
}
return (weightField & Sf);
}
......
......@@ -145,6 +145,21 @@ Note
- take care when using isoSurfaces - these might have duplicate
triangles and so integration might be wrong
Uniformity
\f[
UI(\phi) = 1 - \frac{1}{2 \overline{\phi} A}
\int{\left| W \phi \cdot \hat{n} - \bar{W} \bar{\phi}\right| d\vec{A}}
\,,\;
\bar{\phi} = \frac{\int{W \phi \cdot d\vec{A}}}{\int{W \cdot d\vec{A}}}
\f]
A velocity uniformity index is calculated with no weighting (W=1) and
\f$ \phi = \vec{U} \f$.
A scalar concentration uniformity index is calculated with either
\f$ \rho \vec{U} \f$ or \f$ \vec{U} \f$ for weighting and
\f$ \phi = conc \f$.
See also
Foam::fieldValues
Foam::functionObject
......@@ -232,6 +247,7 @@ public:
opCoV, //!< Coefficient of variation
opAreaNormalAverage, //!< Area average in normal direction
opAreaNormalIntegrate, //!< Area integral in normal direction
opUniformity, //!< Uniformity index
// Weighted variants
......@@ -247,6 +263,9 @@ public:
//! Weighted area integral
opWeightedAreaIntegrate = (opAreaIntegrate | typeWeighted),
//! Weighted uniformity index
opWeightedUniformity = (opUniformity | typeWeighted),
// Variants using absolute weighting
//! Sum using abs weighting
......@@ -261,6 +280,10 @@ public:
//! Area integral using abs weighting
opAbsWeightedAreaIntegrate =
(opWeightedAreaIntegrate | typeAbsolute),
//! Uniformity index using abs weighting
opAbsWeightedUniformity =
(opWeightedUniformity | typeAbsolute),
};
//- Operation type names
......
......@@ -275,10 +275,28 @@ processSameTypeValues
case opAreaNormalAverage:
case opAreaNormalIntegrate:
case opUniformity:
{
// Handled in specializations only
break;
}
case opWeightedUniformity:
case opAbsWeightedUniformity:
{
if (canWeight(weightField))
{
// Change weighting from vector -> scalar and dispatch again
return processValues<Type, scalar>
(
values,
Sf,
weightingFactor(weightField)
);
}
break;
}
}
return result;
......
......@@ -20,23 +20,22 @@ internalField uniform 1000;
boundaryField
{
Default_Boundary_Region
{
type zeroGradient;
}
inlet
{
type fixedValue;
value uniform 1000;
value $internalField;
}
outlet
{
type inletOutlet;
//type zeroGradient;
value uniform 1000;
inletValue uniform 1000;
value $internalField;
inletValue $internalField;;
}
".*"
{
type zeroGradient;
}
}
......
......@@ -33,8 +33,8 @@ boundaryField
outlet
{
type inletOutlet;
value uniform (0 0 0);
inletValue uniform (0 0 0);
value $internalField;
inletValue $internalField;
}
}
......
......@@ -10,7 +10,6 @@ FoamFile
version 2.0;
format ascii;
class volScalarField;
location "0";
object alphat;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -25,17 +24,17 @@ boundaryField
{
type compressible::alphatWallFunction;
Prt 0.85;
value uniform 0;
value $internalField;
}
inlet
{
type calculated;
value uniform 0;
value $internalField;
}
outlet
{
type calculated;
value uniform 0;
value $internalField;
}
}
......
......@@ -10,7 +10,6 @@ FoamFile
version 2.0;
format ascii;
class volScalarField;
location "0";
object epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -27,19 +26,19 @@ boundaryField
Cmu 0.09;
kappa 0.41;
E 9.8;
value uniform 200;
value $internalField;
}
inlet
{
type turbulentMixingLengthDissipationRateInlet;
mixingLength 0.005;
value uniform 200;
value $internalField;
}
outlet
{
type inletOutlet;
inletValue uniform 200;
value uniform 200;
value $internalField;
inletValue $internalField;
}
}
......
......@@ -10,7 +10,6 @@ FoamFile
version 2.0;
format ascii;
class volScalarField;
location "0";
object k;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -24,19 +23,19 @@ boundaryField
Default_Boundary_Region
{
type kqRWallFunction;
value uniform 1;
value $internalField;
}
inlet
{
type turbulentIntensityKineticEnergyInlet;
intensity 0.05;
value uniform 1;
value $internalField;
}
outlet
{
type inletOutlet;
inletValue uniform 1;
value uniform 1;
value $internalField;
inletValue $internalField;
}
}
......
......@@ -10,7 +10,6 @@ FoamFile
version 2.0;
format ascii;
class volScalarField;
location "0";
object nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -27,17 +26,17 @@ boundaryField
Cmu 0.09;
kappa 0.41;
E 9.8;
value uniform 0;
value $internalField;
}
inlet
{
type calculated;
value uniform 0;
value $internalField;
}
outlet
{
type calculated;
value uniform 0;
value $internalField;
}
}
......
......@@ -28,17 +28,17 @@ boundaryField
{
type zeroGradient;
//type mixed;
refValue uniform 110000;
refValue $internalField;
refGradient uniform 0;
valueFraction uniform 0.3;
}
outlet
{
type fixedValue;
value uniform 110000;
value $internalField;
//type mixed;
//refValue uniform 110000;
//refValue $internalField;
//refGradient uniform 0;
//valueFraction uniform 1;
//type transonicOutletPressure;
......@@ -46,7 +46,12 @@ boundaryField
//phi phi;
//gamma 1.4;
//psi psi;
//pInf uniform 110000;
//pInf $internalField;
}
".*"
{
type zeroGradient;
}
}
......
......@@ -47,5 +47,12 @@ graphFormat raw;
runTimeModifiable true;
#include "sampleControls"
functions
{
#include "sampling"
}
// ************************************************************************* //
// -*- C++ -*-
// ************************************************************************* //
// Transcribe volume fields to surfaces.
fieldTransfer
{
type surfMeshes;
libs ("libsampling.so");
log true;
writeControl none;
createOnRead true;
executeControl timeStep;
executeInterval 1;
fields (p rho U T);
derived (rhoU);
_plane
{
type plane;
source cells;
planeType pointAndNormal;
pointAndNormalDict
{
normal (-1 0 0);
point (-0.04 0 0);
}
}
surfaces
(
// Top channel
plane1
{
${_plane}
bounds (-1 0 -1) (0 1 1);
}
// Bottom channel
plane2
{
${_plane}
bounds (-1 -1 -1) (0 0 1);
}
);
}
// ************************************************************************* //
// -*- C++ -*-
// restartTime:
// - a 'one-shot' reset at a particular time
//
// fields [required]
// Pairs of fields to use for calculating the deviation.
// The fields must already exist on the surfaces.
//
// weightField [optional]
// A scalar or vector field for weighting.
//
// postOperation [optional]
// Modify the results by particular operations.
// (none | sqrt)
// The sqrt operation is useful when determining RMS values.
//
// The 'output/write' control triggers the calculation.
__surfaceFieldValue
{
type surfaceFieldValue;
libs ("libfieldFunctionObjects.so");
log on;
enabled true;
writeControl timeStep;
writeInterval 1;
writeFields false;
surfaceFormat vtk;
// writeArea true;
// resetOnStartUp true;
// resetOnOutput false;
// periodicRestart true;
// restartPeriod 0.0005;
}
// ************************************************************************* //
// -*- C++ -*-
// ************************************************************************* //
#include "fieldTransfer"
massflow
{
${__surfaceFieldValue}
regionType surface;
name plane1;
operation areaNormalIntegrate;
fields ( rhoU );
}
areaIntegrate
{
${__surfaceFieldValue}
regionType surface;
name plane1;
operation weightedAreaIntegrate;
weightField rhoU;
fields ( T );
}
// Inflow uniformity
UI1
{
${__surfaceFieldValue}
regionType surface;
name plane1;
operation uniformity;
fields ( U T );
}