diff --git a/etc/bashrc b/etc/bashrc
index 7907ce303a36d6f65bf684701a4be5f121c6180c..0e7c7c62a4b1ff4882b55b14ce339f7ecfc52244 100644
--- a/etc/bashrc
+++ b/etc/bashrc
@@ -137,7 +137,7 @@ _foamSource()
{
while [ $# -ge 1 ]
do
- [ "$FOAM_VERBOSE" -a "$PS1" ] && echo "Sourcing: $1"
+ [ "$FOAM_VERBOSE" -a "$PS1" ] && echo "Sourcing: $1" 1>&2
. $1
shift
done
@@ -155,12 +155,12 @@ _foamEval()
;;
*=)
# name= -> unset name
- [ "$FOAM_VERBOSE" -a "$PS1" ] && echo "unset ${1%=}"
+ [ "$FOAM_VERBOSE" -a "$PS1" ] && echo "unset ${1%=}" 1>&2
eval "unset ${1%=}"
;;
*=*)
# name=value -> export name=value
- [ "$FOAM_VERBOSE" -a "$PS1" ] && echo "export $1"
+ [ "$FOAM_VERBOSE" -a "$PS1" ] && echo "export $1" 1>&2
eval "export $1"
;;
*)
diff --git a/etc/config/aliases.sh b/etc/config/aliases.sh
index 341899e428cb3135aae949ef3209ab515dc6a393..ea5a25a82a6f1735d797ff23756e66da99e74c55 100644
--- a/etc/config/aliases.sh
+++ b/etc/config/aliases.sh
@@ -53,7 +53,7 @@ unset foamPV
foamPV()
{
. $WM_PROJECT_DIR/etc/config/paraview.sh ParaView_VERSION=$1
- echo "paraview-$ParaView_VERSION (major: $ParaView_MAJOR)"
+ echo "paraview-$ParaView_VERSION (major: $ParaView_MAJOR)" 1>&2
}
diff --git a/etc/config/settings.sh b/etc/config/settings.sh
index 1ece2c9cff34fe8958df364a488fef09ae5b3087..1a7e8237f9da1b45abb08e70adc8ccab3ac13852 100644
--- a/etc/config/settings.sh
+++ b/etc/config/settings.sh
@@ -93,7 +93,7 @@ Linux)
export WM_LDFLAGS='-m64'
;;
*)
- echo "Unknown WM_ARCH_OPTION '$WM_ARCH_OPTION', should be 32 or 64"
+ echo "Unknown WM_ARCH_OPTION '$WM_ARCH_OPTION', should be 32 or 64" 1>&2
;;
esac
;;
@@ -135,7 +135,7 @@ Linux)
;;
*)
- echo Unknown processor type `uname -m` for Linux
+ echo Unknown processor type `uname -m` for Linux 1>&2
;;
esac
;;
@@ -152,7 +152,7 @@ SunOS)
;;
*) # an unsupported operating system
- cat <<USAGE
+ /bin/cat <<USAGE 1>&2
Your "$WM_ARCH" operating system is not supported by this release
of OpenFOAM. For further assistance, please contact www.OpenFOAM.org
@@ -268,11 +268,11 @@ OpenFOAM | ThirdParty)
#clang_version=llvm-svn
;;
*)
- echo
- echo "Warning in $WM_PROJECT_DIR/etc/config/settings.sh:"
- echo " Unknown OpenFOAM compiler type '$WM_COMPILER'"
- echo " Please check your settings"
- echo
+ echo 1>&2
+ echo "Warning in $WM_PROJECT_DIR/etc/config/settings.sh:" 1>&2
+ echo " Unknown OpenFOAM compiler type '$WM_COMPILER'" 1>&2
+ echo " Please check your settings" 1>&2
+ echo 1>&2
;;
esac
@@ -288,11 +288,11 @@ OpenFOAM | ThirdParty)
# Check that the compiler directory can be found
[ -d "$gccDir" ] || {
- echo
- echo "Warning in $WM_PROJECT_DIR/etc/config/settings.sh:"
- echo " Cannot find $gccDir installation."
- echo " Please install this compiler version or if you wish to use the system compiler,"
- echo " change the 'foamCompiler' setting to 'system'"
+ echo 1>&2
+ echo "Warning in $WM_PROJECT_DIR/etc/config/settings.sh:" 1>&2
+ echo " Cannot find $gccDir installation." 1>&2
+ echo " Please install this compiler version or if you wish to use the system compiler," 1>&2
+ echo " change the 'foamCompiler' setting to 'system'" 1>&2
echo
}
@@ -325,12 +325,12 @@ OpenFOAM | ThirdParty)
# Check that the compiler directory can be found
[ -d "$clangDir" ] || {
- echo
- echo "Warning in $WM_PROJECT_DIR/etc/config/settings.sh:"
- echo " Cannot find $clangDir installation."
- echo " Please install this compiler version or if you wish to use the system compiler,"
- echo " change the 'foamCompiler' setting to 'system'"
- echo
+ echo 1>&2
+ echo "Warning in $WM_PROJECT_DIR/etc/config/settings.sh:" 1>&2
+ echo " Cannot find $clangDir installation." 1>&2
+ echo " Please install this compiler version or if you wish to use the system compiler," 1>&2
+ echo " change the 'foamCompiler' setting to 'system'" 1>&2
+ echo 1>&2
}
_foamAddMan $clangDir/share/man
@@ -415,10 +415,10 @@ SYSTEMOPENMPI)
if [ "$FOAM_VERBOSE" -a "$PS1" ]
then
- echo "Using system installed MPI:"
- echo " compile flags : $PINC"
- echo " link flags : $PLIBS"
- echo " libmpi dir : $libDir"
+ echo "Using system installed MPI:" 1>&2
+ echo " compile flags : $PINC" 1>&2
+ echo " link flags : $PLIBS" 1>&2
+ echo " libmpi dir : $libDir" 1>&2
fi
_foamAddLib $libDir
@@ -492,7 +492,7 @@ HPMPI)
_foamAddLib $MPI_ARCH_PATH/lib/linux_ia64
;;
*)
- echo Unknown processor type `uname -m` for Linux
+ echo Unknown processor type `uname -m` 1>&2
;;
esac
;;
@@ -526,13 +526,10 @@ QSMPI)
;;
SGIMPI)
- lastCharID=$(( ${#MPI_ROOT} - 1 ))
- if [ "${MPI_ROOT:$lastCharID:1}" == '/' ]
- then
- MPI_ROOT=${MPI_ROOT:0:$lastCharID}
- fi
+ # no trailing slash
+ [ "${MPI_ROOT%/}" = "${MPI_ROOT}" ] || MPI_ROOT="${MPI_ROOT%/}"
- export FOAM_MPI=${MPI_ROOT##*/}
+ export FOAM_MPI="${MPI_ROOT##*/}"
export MPI_ARCH_PATH=$MPI_ROOT
if [ ! -d "$MPI_ROOT" -o -z "$MPI_ARCH_PATH" ]
@@ -545,9 +542,9 @@ SGIMPI)
if [ "$FOAM_VERBOSE" -a "$PS1" ]
then
- echo "Using SGI MPT:"
- echo " MPI_ROOT : $MPI_ROOT"
- echo " FOAM_MPI : $FOAM_MPI"
+ echo "Using SGI MPT:" 1>&2
+ echo " MPI_ROOT : $MPI_ROOT" 1>&2
+ echo " FOAM_MPI : $FOAM_MPI" 1>&2
fi
_foamAddPath $MPI_ARCH_PATH/bin
@@ -555,13 +552,10 @@ SGIMPI)
;;
INTELMPI)
- lastCharID=$(( ${#MPI_ROOT} - 1 ))
- if [ "${MPI_ROOT:$lastCharID:1}" == '/' ]
- then
- MPI_ROOT=${MPI_ROOT:0:$lastCharID}
- fi
+ # no trailing slash
+ [ "${MPI_ROOT%/}" = "${MPI_ROOT}" ] || MPI_ROOT="${MPI_ROOT%/}"
- export FOAM_MPI=${MPI_ROOT##*/}
+ export FOAM_MPI="${MPI_ROOT##*/}"
export MPI_ARCH_PATH=$MPI_ROOT
if [ ! -d "$MPI_ROOT" -o -z "$MPI_ARCH_PATH" ]
@@ -574,9 +568,9 @@ INTELMPI)
if [ "$FOAM_VERBOSE" -a "$PS1" ]
then
- echo "Using INTEL MPI:"
- echo " MPI_ROOT : $MPI_ROOT"
- echo " FOAM_MPI : $FOAM_MPI"
+ echo "Using INTEL MPI:" 1>&2
+ echo " MPI_ROOT : $MPI_ROOT" 1>&2
+ echo " FOAM_MPI : $FOAM_MPI" 1>&2
fi
_foamAddPath $MPI_ARCH_PATH/bin64
diff --git a/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.C b/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.C
index c3f042fe0160277a4773b46da6f4b334259e2e41..2b73a471f48ad3fbf57a8dd97645ae1cc6ba51f6 100644
--- a/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.C
+++ b/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.C
@@ -369,7 +369,7 @@ void Foam::rotorDiskSource::constructGeometry()
// blade flap angle [radians]
scalar beta =
- flap_.beta0 - flap_.beta1*cos(psi) - flap_.beta2*sin(psi);
+ flap_.beta0 - flap_.beta1c*cos(psi) - flap_.beta2s*sin(psi);
// determine rotation tensor to convert from planar system into the
// rotor cone system
@@ -659,11 +659,11 @@ bool Foam::rotorDiskSource::read(const dictionary& dict)
const dictionary& flapCoeffs(coeffs_.subDict("flapCoeffs"));
flapCoeffs.lookup("beta0") >> flap_.beta0;
- flapCoeffs.lookup("beta1") >> flap_.beta1;
- flapCoeffs.lookup("beta2") >> flap_.beta2;
+ flapCoeffs.lookup("beta1c") >> flap_.beta1c;
+ flapCoeffs.lookup("beta2s") >> flap_.beta2s;
flap_.beta0 = degToRad(flap_.beta0);
- flap_.beta1 = degToRad(flap_.beta1);
- flap_.beta2 = degToRad(flap_.beta2);
+ flap_.beta1c = degToRad(flap_.beta1c);
+ flap_.beta2s = degToRad(flap_.beta2s);
// create co-ordinate system
diff --git a/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.H b/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.H
index 7fd7d33964602b8574da194ff178eb882b546445..17deb44816be960353e83e5b8ee118ab68adb92d 100644
--- a/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.H
+++ b/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.H
@@ -51,8 +51,8 @@ Description
flapCoeffs
{
beta0 0; // coning angle [deg]
- beta1 0; // lateral flapping coeff
- beta2 0; // longitudinal flapping coeff
+ beta1c 0; // lateral flapping coeff (cos coeff)
+ beta2s 0; // longitudinal flapping coeff (sin coeff)
}
blade
@@ -78,7 +78,6 @@ Description
local : use local flow conditions
-
SourceFiles
rotorDiskSource.C
rotorDiskSourceTemplates.C
@@ -138,8 +137,8 @@ protected:
struct flapData
{
scalar beta0; // coning angle
- scalar beta1; // lateral flapping coeff
- scalar beta2; // longitudinal flapping coeff
+ scalar beta1c; // lateral flapping coeff (cos coeff)
+ scalar beta2s; // longitudinal flapping coeff (sin coeff)
};
@@ -148,7 +147,7 @@ protected:
//- Name of density field
word rhoName_;
- //- Reference density if rhoName = 'none'
+ //- Reference density for rhoName = 'none'
scalar rhoRef_;
//- Rotational speed [rad/s]
@@ -254,6 +253,9 @@ public:
// Access
+ //- Return the reference density for rhoName = 'none'
+ inline scalar rhoRef() const;
+
//- Return the rotational speed [rad/s]
// Positive anti-clockwise when looking along -ve lift direction
inline scalar omega() const;
diff --git a/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSourceI.H b/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSourceI.H
index 7c194e8713e4beb538292f9e2fd3a6f1c44a9ed5..9d095e846b9d3343bcef811a3301de304ac034fa 100644
--- a/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSourceI.H
+++ b/src/fieldSources/basicSource/rotorDiskSource/rotorDiskSourceI.H
@@ -27,11 +27,18 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+Foam::scalar Foam::rotorDiskSource::rhoRef() const
+{
+ return rhoRef_;
+}
+
+
Foam::scalar Foam::rotorDiskSource::omega() const
{
return omega_;
}
+
const Foam::List<Foam::point>& Foam::rotorDiskSource::x() const
{
return x_;
diff --git a/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.C b/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.C
index d53d73f3ac454741e6c81556f8fc29bbdcb97df4..4ed3f7712019026712f3ed91c07e5dc42d0e8605 100644
--- a/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.C
+++ b/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.C
@@ -70,25 +70,33 @@ Foam::vector Foam::targetCoeffTrim::calcCoeffs
{
label cellI = cells[i];
- // create normalisation coefficient
- scalar radius = x[i].x();
- scalar coeff2 = coeff1*pow4(radius);
- if (compressible)
- {
- coeff2 *= trho()[cellI];
- }
-
vector fc = force[cellI];
vector mc = fc^(C[cellI] - origin);
- // add to coefficient vector
- cf[0] += (fc & yawAxis)/(coeff2 + ROOTVSMALL);
- cf[1] += (mc & pitchAxis)/(coeff2*radius + ROOTVSMALL);
- cf[2] += (mc & rollAxis)/(coeff2*radius + ROOTVSMALL);
+ if (useCoeffs_)
+ {
+ scalar radius = x[i].x();
+ scalar coeff2 = coeff1*pow4(radius);
+ if (compressible)
+ {
+ coeff2 *= trho()[cellI];
+ }
+
+ // add to coefficient vector
+ cf[0] += (fc & yawAxis)/(coeff2 + ROOTVSMALL);
+ cf[1] += (mc & pitchAxis)/(coeff2*radius + ROOTVSMALL);
+ cf[2] += (mc & rollAxis)/(coeff2*radius + ROOTVSMALL);
+ }
+ else
+ {
+ cf[0] += fc & yawAxis;
+ cf[1] += mc & pitchAxis;
+ cf[2] += mc & rollAxis;
+ }
}
reduce(cf, sumOp<vector>());
-
+
return cf;
}
@@ -103,6 +111,7 @@ Foam::targetCoeffTrim::targetCoeffTrim
:
trimModel(rotor, dict, typeName),
calcFrequency_(-1),
+ useCoeffs_(true),
target_(vector::zero),
theta_(vector::zero),
nIter_(50),
@@ -128,9 +137,16 @@ void Foam::targetCoeffTrim::read(const dictionary& dict)
trimModel::read(dict);
const dictionary& targetDict(coeffs_.subDict("target"));
- target_[0] = readScalar(targetDict.lookup("thrustCoeff"));
- target_[1] = readScalar(targetDict.lookup("pitchCoeff"));
- target_[2] = readScalar(targetDict.lookup("rollCoeff"));
+ useCoeffs_ = targetDict.lookupOrDefault<bool>("useCoeffs", true);
+ word ext = "";
+ if (useCoeffs_)
+ {
+ ext = "Coeff";
+ }
+
+ target_[0] = readScalar(targetDict.lookup("thrust" + ext));
+ target_[1] = readScalar(targetDict.lookup("pitch" + ext));
+ target_[2] = readScalar(targetDict.lookup("roll" + ext));
const dictionary& pitchAngleDict(coeffs_.subDict("pitchAngles"));
theta_[0] = degToRad(readScalar(pitchAngleDict.lookup("theta0Ini")));
@@ -173,8 +189,16 @@ void Foam::targetCoeffTrim::correct(const vectorField& U, vectorField& force)
{
if (rotor_.mesh().time().timeIndex() % calcFrequency_ == 0)
{
+ word calcType = "forces";
+ if (useCoeffs_)
+ {
+ calcType = "coefficients";
+ }
+
Info<< type() << ":" << nl
- << " solving for target trim coefficients" << nl;
+ << " solving for target trim " << calcType << nl;
+
+ const scalar rhoRef = rotor_.rhoRef();
// iterate to find new pitch angles to achieve target force
scalar err = GREAT;
@@ -210,7 +234,7 @@ void Foam::targetCoeffTrim::correct(const vectorField& U, vectorField& force)
}
// calculate the change in pitch angle vector
- vector dt = inv(J) & (target_ - old);
+ vector dt = inv(J) & (target_/rhoRef - old);
// update pitch angles
vector thetaNew = theta_ + relax_*dt;
@@ -235,10 +259,10 @@ void Foam::targetCoeffTrim::correct(const vectorField& U, vectorField& force)
<< "), iterations = " << iter << endl;
}
- Info<< " current and target coefficients:" << nl
- << " thrust = " << old[0] << ", " << target_[0] << nl
- << " pitch = " << old[1] << ", " << target_[1] << nl
- << " roll = " << old[2] << ", " << target_[2] << nl
+ Info<< " current and target " << calcType << nl
+ << " thrust = " << old[0]*rhoRef << ", " << target_[0] << nl
+ << " pitch = " << old[1]*rhoRef << ", " << target_[1] << nl
+ << " roll = " << old[2]*rhoRef << ", " << target_[2] << nl
<< " new pitch angles [deg]:" << nl
<< " theta0 = " << radToDeg(theta_[0]) << nl
<< " theta1c = " << radToDeg(theta_[1]) << nl
diff --git a/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.H b/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.H
index e387a4c1630fa36b14cdd80e76d59eb45bb28391..24d5732e99f719c5053f0dbdac6063259d665022 100644
--- a/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.H
+++ b/src/fieldSources/basicSource/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.H
@@ -25,7 +25,7 @@ Class
Foam::targetCoeffTrim
Description
- Target trim coefficients
+ Target trim forces/coefficients
Solves:
@@ -43,8 +43,8 @@ Description
Newton-Raphson iterative method. The solver tolerance can be user-input,
using the 'tol' entry.
-
- The coefficients are determined using:
+ If coefficients are requested (useCoeffs = true), the force and moments
+ are normalised using:
force
c = ---------------------------------
@@ -97,6 +97,9 @@ protected:
//- Number of iterations between calls to 'correct'
label calcFrequency_;
+ //- Flag to indicate whether to solve coeffs (true) or forces (false)
+ bool useCoeffs_;
+
//- Target coefficient vector (thrust force, roll moment, pitch moment)
vector target_;
diff --git a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
index 8e291b75b7185d431304237a3040e6feed51320c..3381bbc087f2f82e036bc72b503a2b59e696eb78 100644
--- a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
+++ b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.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-2012 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -246,7 +246,7 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
const scalar mass = p.mass();
const forceSuSp Fcp = forces.calcCoupled(p, dt, mass, Re, muAv);
const forceSuSp Fncp = forces.calcNonCoupled(p, dt, mass, Re, muAv);
- scalar tMom = 1.0/(Fcp.Sp() + Fncp.Sp());
+ scalar tMom = mass/(Fcp.Sp() + Fncp.Sp());
const vector g = td.cloud().g().value();
@@ -301,7 +301,7 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
child->tc() = -GREAT;
child->ms() = 0.0;
child->injector() = this->injector();
- child->tMom() = 1.0/(Fcp.Sp() + Fncp.Sp());
+ child->tMom() = massChild/(Fcp.Sp() + Fncp.Sp());
child->user() = 0.0;
child->setCellValues(td, dt, cellI);
diff --git a/src/postProcessing/functionObjects/field/Make/files b/src/postProcessing/functionObjects/field/Make/files
index e351784bdcb4e26fa40de04f36d6995f51f9ccf3..0abf314f96739f68133c1013b8f81ead28c17116 100644
--- a/src/postProcessing/functionObjects/field/Make/files
+++ b/src/postProcessing/functionObjects/field/Make/files
@@ -41,4 +41,7 @@ wallBoundedStreamLine/wallBoundedParticle.C
surfaceInterpolateFields/surfaceInterpolateFields.C
surfaceInterpolateFields/surfaceInterpolateFieldsFunctionObject.C
+regionSizeDistribution/regionSizeDistribution.C
+regionSizeDistribution/regionSizeDistributionFunctionObject.C
+
LIB = $(FOAM_LIBBIN)/libfieldFunctionObjects
diff --git a/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistribution.C b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistribution.C
new file mode 100644
index 0000000000000000000000000000000000000000..b3aaac0dbf9fa2b525d96654f0fd494a470a449b
--- /dev/null
+++ b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistribution.C
@@ -0,0 +1,588 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "regionSizeDistribution.H"
+#include "volFields.H"
+#include "regionSplit.H"
+#include "fvcVolumeIntegrate.H"
+#include "Histogram.H"
+#include "mathematicalConstants.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(regionSizeDistribution, 0);
+
+ //- plus op for FixedList<scalar>
+ template<class T, unsigned Size>
+ class ListPlusEqOp
+ {
+ public:
+ void operator()
+ (
+ FixedList<T, Size>& x,
+ const FixedList<T, Size>& y
+ ) const
+ {
+ forAll(x, i)
+ {
+ x[i] += y[i];
+ }
+ }
+ };
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::regionSizeDistribution::regionSizeDistribution
+(
+ const word& name,
+ const objectRegistry& obr,
+ const dictionary& dict,
+ const bool loadFromFiles
+)
+:
+ name_(name),
+ obr_(obr),
+ active_(true),
+ alphaName_(dict.lookup("field")),
+ patchNames_(dict.lookup("patches"))
+{
+ // Check if the available mesh is an fvMesh, otherwise deactivate
+ if (!isA<fvMesh>(obr_))
+ {
+ active_ = false;
+ WarningIn
+ (
+ "regionSizeDistribution::regionSizeDistribution"
+ "(const objectRegistry&, const dictionary&)"
+ ) << "No fvMesh available, deactivating." << nl
+ << endl;
+ }
+
+ read(dict);
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::regionSizeDistribution::~regionSizeDistribution()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::regionSizeDistribution::read(const dictionary& dict)
+{
+ if (active_)
+ {
+ dict.lookup("field") >> alphaName_;
+ dict.lookup("patches") >> patchNames_;
+ dict.lookup("threshold") >> threshold_;
+ dict.lookup("volFraction") >> volFraction_;
+ dict.lookup("nBins") >> nBins_;
+
+ word format(dict.lookup("setFormat"));
+ formatterPtr_ = writer<scalar>::New(format);
+ }
+}
+
+
+void Foam::regionSizeDistribution::execute()
+{
+ // Do nothing - only valid on write
+}
+
+
+void Foam::regionSizeDistribution::end()
+{
+ // Do nothing - only valid on write
+}
+
+
+void Foam::regionSizeDistribution::write()
+{
+ if (active_)
+ {
+ const fvMesh& mesh = refCast<const fvMesh>(obr_);
+
+ autoPtr<volScalarField> alphaPtr;
+ if (obr_.foundObject<volScalarField>(alphaName_))
+ {
+ Info<< "Looking up field " << alphaName_ << endl;
+ }
+ else
+ {
+ Info<< "Reading field " << alphaName_ << endl;
+ alphaPtr.reset
+ (
+ new volScalarField
+ (
+ IOobject
+ (
+ alphaName_,
+ mesh.time().timeName(),
+ mesh,
+ IOobject::MUST_READ,
+ IOobject::NO_WRITE
+ ),
+ mesh
+ )
+ );
+ }
+
+
+ const volScalarField& alpha =
+ (
+ alphaPtr.valid()
+ ? alphaPtr()
+ : obr_.lookupObject<volScalarField>(alphaName_)
+ );
+
+ Info<< "Volume of alpha = "
+ << fvc::domainIntegrate(alpha).value()
+ << endl;
+
+ const scalar meshVol = gSum(mesh.V());
+ Info<< "Mesh volume = " << meshVol << endl;
+ Info<< "Background region volume limit = " << volFraction_*meshVol
+ << endl;
+
+
+ // Determine blocked faces
+ boolList blockedFace(mesh.nFaces(), false);
+ label nBlocked = 0;
+
+ {
+ for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
+ {
+ scalar ownVal = alpha[mesh.faceOwner()[faceI]];
+ scalar neiVal = alpha[mesh.faceNeighbour()[faceI]];
+
+ if
+ (
+ (ownVal < threshold_ && neiVal > threshold_)
+ || (ownVal > threshold_ && neiVal < threshold_)
+ )
+ {
+ blockedFace[faceI] = true;
+ nBlocked++;
+ }
+ }
+
+ // Block coupled faces
+ forAll(alpha.boundaryField(), patchI)
+ {
+ const fvPatchScalarField& fvp = alpha.boundaryField()[patchI];
+ if (fvp.coupled())
+ {
+ tmp<scalarField> townFld(fvp.patchInternalField());
+ const scalarField& ownFld = townFld();
+ tmp<scalarField> tnbrFld(fvp.patchNeighbourField());
+ const scalarField& nbrFld = tnbrFld();
+
+ label start = fvp.patch().patch().start();
+
+ forAll(ownFld, i)
+ {
+ scalar ownVal = ownFld[i];
+ scalar neiVal = nbrFld[i];
+
+ if
+ (
+ (ownVal < threshold_ && neiVal > threshold_)
+ || (ownVal > threshold_ && neiVal < threshold_)
+ )
+ {
+ blockedFace[start+i] = true;
+ nBlocked++;
+ }
+ }
+ }
+ }
+ }
+
+
+ regionSplit regions(mesh, blockedFace);
+
+ Info<< "Determined " << regions.nRegions() << " disconnected regions"
+ << endl;
+
+
+ if (debug)
+ {
+ volScalarField region
+ (
+ IOobject
+ (
+ "region",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::NO_READ,
+ IOobject::NO_WRITE
+ ),
+ mesh,
+ dimensionedScalar("zero", dimless, 0)
+ );
+ Info<< "Dumping region as volScalarField to " << region.name()
+ << endl;
+
+ forAll(regions, cellI)
+ {
+ region[cellI] = regions[cellI];
+ }
+ region.correctBoundaryConditions();
+ region.write();
+ }
+
+
+ // Sum all regions
+ Map<Pair<scalar> > regionVolume(regions.nRegions()/Pstream::nProcs());
+ forAll(alpha, cellI)
+ {
+ scalar cellVol = mesh.V()[cellI];
+ scalar alphaVol = alpha[cellI]*cellVol;
+
+ label regionI = regions[cellI];
+
+ Map<Pair<scalar> >::iterator fnd = regionVolume.find(regionI);
+ if (fnd == regionVolume.end())
+ {
+ regionVolume.insert
+ (
+ regionI,
+ Pair<scalar>(cellVol, alphaVol)
+ );
+ }
+ else
+ {
+ fnd().first() += cellVol;
+ fnd().second() += alphaVol;
+ }
+ }
+ Pstream::mapCombineGather(regionVolume, ListPlusEqOp<scalar, 2>());
+ Pstream::mapCombineScatter(regionVolume);
+
+
+ if (debug)
+ {
+ Info<< token::TAB << "Region"
+ << token::TAB << "Volume(mesh)"
+ << token::TAB << "Volume(" << alpha.name() << "):"
+ << endl;
+ scalar meshSumVol = 0.0;
+ scalar alphaSumVol = 0.0;
+
+ forAllConstIter(Map<Pair<scalar> >, regionVolume, iter)
+ {
+ Info<< token::TAB << iter.key()
+ << token::TAB << iter().first()
+ << token::TAB << iter().second() << endl;
+
+ meshSumVol += iter().first();
+ alphaSumVol += iter().second();
+ }
+ Info<< token::TAB << "Total:"
+ << token::TAB << meshSumVol
+ << token::TAB << alphaSumVol << endl;
+ Info<< endl;
+ }
+
+
+
+ // Mark all regions starting at patches
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // Count number of patch faces (just for initial sizing)
+ label nPatchFaces = 0;
+ forAll(patchNames_, i)
+ {
+ const word& pName = patchNames_[i];
+ label patchI = mesh.boundaryMesh().findPatchID(pName);
+ if (patchI == -1)
+ {
+ WarningIn("regionSizeDistribution::write()")
+ << "Cannot find patch " << pName << ". Valid patches are "
+ << mesh.boundaryMesh().names()
+ << endl;
+ }
+ else
+ {
+ nPatchFaces += mesh.boundaryMesh()[patchI].size();
+ }
+ }
+
+ Map<label> keepRegions(nPatchFaces);
+ forAll(patchNames_, i)
+ {
+ const word& pName = patchNames_[i];
+
+ label patchI = mesh.boundaryMesh().findPatchID(pName);
+ if (patchI != -1)
+ {
+ const polyPatch& pp = mesh.boundaryMesh()[patchI];
+
+ // Collect all regions on the patch
+ const labelList& faceCells = pp.faceCells();
+
+ forAll(faceCells, i)
+ {
+ keepRegions.insert
+ (
+ regions[faceCells[i]],
+ Pstream::myProcNo()
+ );
+ }
+ }
+ }
+
+
+ // Make sure all the processors have the same set of regions
+ Pstream::mapCombineGather(keepRegions, minEqOp<label>());
+ Pstream::mapCombineScatter(keepRegions);
+
+ Info<< "Patch connected regions (liquid core):" << endl;
+ forAllConstIter(Map<label>, keepRegions, iter)
+ {
+ label regionI = iter.key();
+ Pair<scalar>& vols = regionVolume[regionI];
+ Info<< token::TAB << iter.key()
+ << token::TAB << vols.first()
+ << token::TAB << vols.second() << endl;
+
+ }
+ Info<< endl;
+
+ Info<< "Background regions:" << endl;
+ forAllConstIter(Map<Pair<scalar> >, regionVolume, iter)
+ {
+ if
+ (
+ !keepRegions.found(iter.key())
+ && iter().first() >= volFraction_*meshVol
+ )
+ {
+ Info<< token::TAB << iter.key()
+ << token::TAB << iter().first()
+ << token::TAB << iter().second() << endl;
+ }
+ }
+ Info<< endl;
+
+
+ // Split alpha field
+ // ~~~~~~~~~~~~~~~~~
+ // Split into
+ // - liquidCore : region connected to inlet patches
+ // - per region a volume : for all other regions
+ // - backgroundAlpha : remaining alpha
+
+
+ // Construct field
+ volScalarField liquidCore
+ (
+ IOobject
+ (
+ alphaName_ + "_liquidCore",
+ obr_.time().timeName(),
+ obr_,
+ IOobject::NO_READ
+ ),
+ alpha,
+ fvPatchField<scalar>::calculatedType()
+ );
+
+ volScalarField backgroundAlpha
+ (
+ IOobject
+ (
+ alphaName_ + "_background",
+ obr_.time().timeName(),
+ obr_,
+ IOobject::NO_READ
+ ),
+ alpha,
+ fvPatchField<scalar>::calculatedType()
+ );
+
+
+ // Knock out any cell not in keepRegions
+ forAll(liquidCore, cellI)
+ {
+ label regionI = regions[cellI];
+ if (keepRegions.found(regionI))
+ {
+ backgroundAlpha[cellI] = 0;
+ }
+ else
+ {
+ liquidCore[cellI] = 0;
+
+ scalar regionVol = regionVolume[regionI].first();
+ if (regionVol < volFraction_*meshVol)
+ {
+ backgroundAlpha[cellI] = 0;
+ }
+ }
+ }
+ liquidCore.correctBoundaryConditions();
+ backgroundAlpha.correctBoundaryConditions();
+
+ Info<< "Volume of liquid-core = "
+ << fvc::domainIntegrate(liquidCore).value()
+ << endl;
+
+ Info<< "Writing liquid-core field to " << liquidCore.name() << endl;
+ liquidCore.write();
+
+ Info<< "Volume of background = "
+ << fvc::domainIntegrate(backgroundAlpha).value()
+ << endl;
+
+ Info<< "Writing background field to " << backgroundAlpha.name() << endl;
+ backgroundAlpha.write();
+
+
+
+ // Collect histogram
+ if (Pstream::master())
+ {
+ DynamicList<scalar> diameters(regionVolume.size());
+ forAllConstIter(Map<Pair<scalar> >, regionVolume, iter)
+ {
+ if (!keepRegions.found(iter.key()))
+ {
+ if (iter().first() < volFraction_*meshVol)
+ {
+ scalar v = iter().second();
+ //scalar diam = Foam::cbrt(v*6/mathematicalConstant::pi);
+ scalar diam =
+ Foam::cbrt(v*6/constant::mathematical::pi);
+ diameters.append(diam);
+ }
+ }
+ }
+
+ if (diameters.size())
+ {
+ scalar maxDiam = max(diameters);
+ scalar minDiam = 0.0;
+
+ Info<< "Maximum diameter:" << maxDiam << endl;
+
+ Histogram<List<scalar> > bins
+ (
+ minDiam,
+ maxDiam,
+ nBins_,
+ diameters
+ );
+
+ /* 1.7.x
+ scalarField xBin(nBins_);
+
+ scalar dx = (maxDiam-minDiam)/nBins_;
+ scalar x = 0.5*dx;
+ forAll(bins.counts(), i)
+ {
+ xBin[i] = x;
+ x += dx;
+ }
+
+ scalarField normalisedCount(bins.counts().size());
+ forAll(bins.counts(), i)
+ {
+ normalisedCount[i] = 1.0*bins.counts()[i];
+ }
+
+ const coordSet coords
+ (
+ "diameter",
+ "x",
+ xBin
+ );
+ */
+
+ pointField xBin(nBins_);
+ scalar dx = (maxDiam - minDiam)/nBins_;
+ scalar x = 0.5*dx;
+ forAll(bins.counts(), i)
+ {
+ xBin[i] = point(x, 0, 0);
+ x += dx;
+ }
+
+ scalarField normalisedCount(bins.counts().size());
+ forAll(bins.counts(), i)
+ {
+ normalisedCount[i] = 1.0*bins.counts()[i];
+ }
+
+ const coordSet coords
+ (
+ "diameter",
+ "x",
+ xBin,
+ mag(xBin)
+ );
+ const wordList valNames(1, "count");
+
+
+ fileName outputPath;
+ if (Pstream::parRun())
+ {
+ outputPath = mesh.time().path()/".."/name_;
+ }
+ else
+ {
+ outputPath = mesh.time().path()/name_;
+ }
+
+ if (mesh.name() != fvMesh::defaultRegion)
+ {
+ outputPath = outputPath/mesh.name();
+ }
+
+ mkDir(outputPath/mesh.time().timeName());
+ OFstream str
+ (
+ outputPath
+ / mesh.time().timeName()
+ / formatterPtr_().getFileName(coords, valNames)
+ );
+ Info<< "Writing distribution to " << str.name() << endl;
+
+ List<const scalarField*> valPtrs(1);
+ valPtrs[0] = &normalisedCount;
+ formatterPtr_().write(coords, valNames, valPtrs, str);
+ }
+ }
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistribution.H b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistribution.H
new file mode 100644
index 0000000000000000000000000000000000000000..033b5f9aeca13c2e81f26f3b0fe52d32b4b90b7d
--- /dev/null
+++ b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistribution.H
@@ -0,0 +1,161 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::regionSizeDistribution
+
+Description
+ Looks up a field, interpolates it to the faces and determines a connected
+ region from a patch where the field is above a certain value.
+ - Writes a field containing all regions starting at given patch
+ ('liquid core')
+ - All other regions are summed for volume and a histogram is calculated.
+
+SourceFiles
+ regionSizeDistribution.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef regionSizeDistribution_H
+#define regionSizeDistribution_H
+
+#include "pointFieldFwd.H"
+#include "writer.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward declaration of classes
+class objectRegistry;
+class dictionary;
+class mapPolyMesh;
+
+/*---------------------------------------------------------------------------*\
+ Class regionSizeDistribution Declaration
+\*---------------------------------------------------------------------------*/
+
+class regionSizeDistribution
+{
+ // Private data
+
+ //- Name of this set of regionSizeDistribution objects
+ word name_;
+
+ const objectRegistry& obr_;
+
+ //- on/off switch
+ bool active_;
+
+ //- Name of field
+ word alphaName_;
+
+ //- Patches to walk from
+ wordList patchNames_;
+
+ //- Clip value
+ scalar threshold_;
+
+ //- Background region volFraction
+ scalar volFraction_;
+
+ //- Mumber of bins
+ label nBins_;
+
+ //- Output formatter to write
+ autoPtr<writer<scalar> > formatterPtr_;
+
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ regionSizeDistribution(const regionSizeDistribution&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const regionSizeDistribution&);
+
+
+public:
+
+ //- Runtime type information
+ TypeName("regionSizeDistribution");
+
+
+ // Constructors
+
+ //- Construct for given objectRegistry and dictionary.
+ // Allow the possibility to load fields from files
+ regionSizeDistribution
+ (
+ const word& name,
+ const objectRegistry&,
+ const dictionary&,
+ const bool loadFromFiles = false
+ );
+
+
+ // Destructor
+
+ virtual ~regionSizeDistribution();
+
+
+ // Member Functions
+
+ //- Return name of the set of regionSizeDistribution
+ virtual const word& name() const
+ {
+ return name_;
+ }
+
+ //- Read the regionSizeDistribution data
+ virtual void read(const dictionary&);
+
+ //- Execute, currently does nothing
+ virtual void execute();
+
+ //- Execute at the final time-loop, currently does nothing
+ virtual void end();
+
+ //- Calculate the regionSizeDistribution and write
+ virtual void write();
+
+ //- Update for changes of mesh
+ virtual void updateMesh(const mapPolyMesh&)
+ {}
+
+ //- Update for changes of mesh
+ virtual void movePoints(const pointField&)
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistributionFunctionObject.C b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistributionFunctionObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..21e25fd6b2d85ea0e221672857fb5892bb6ea9e5
--- /dev/null
+++ b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistributionFunctionObject.C
@@ -0,0 +1,46 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "regionSizeDistributionFunctionObject.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineNamedTemplateTypeNameAndDebug
+ (
+ regionSizeDistributionFunctionObject,
+ 0
+ );
+
+ addToRunTimeSelectionTable
+ (
+ functionObject,
+ regionSizeDistributionFunctionObject,
+ dictionary
+ );
+}
+
+// ************************************************************************* //
diff --git a/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistributionFunctionObject.H b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistributionFunctionObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..97526d0c3ab3e71dd69941faa27786ffbc4b3b8b
--- /dev/null
+++ b/src/postProcessing/functionObjects/field/regionSizeDistribution/regionSizeDistributionFunctionObject.H
@@ -0,0 +1,54 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Typedef
+ Foam::regionSizeDistributionFunctionObject
+
+Description
+ FunctionObject wrapper around regionSizeDistribution to allow it to be
+ created via the functions list within controlDict.
+
+SourceFiles
+ regionSizeDistributionFunctionObject.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef regionSizeDistributionFunctionObject_H
+#define regionSizeDistributionFunctionObject_H
+
+#include "regionSizeDistribution.H"
+#include "OutputFilterFunctionObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ typedef OutputFilterFunctionObject<regionSizeDistribution>
+ regionSizeDistributionFunctionObject;
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/turbulenceModels/incompressible/LES/Make/files b/src/turbulenceModels/incompressible/LES/Make/files
index ff2bde02f26acb7c6c843365b92b3883df737b3f..eadcf1a99f7aee73632187963de14f5c6470fbe9 100644
--- a/src/turbulenceModels/incompressible/LES/Make/files
+++ b/src/turbulenceModels/incompressible/LES/Make/files
@@ -34,5 +34,4 @@ wallFunctions=derivedFvPatchFields/wallFunctions
nuSgsWallFunctions=$(wallFunctions)/nuSgsWallFunctions
$(nuSgsWallFunctions)/nuSgsUSpaldingWallFunction/nuSgsUSpaldingWallFunctionFvPatchScalarField.C
-
LIB = $(FOAM_LIBBIN)/libincompressibleLESModels