From a3788fe854059e0303b90bae9cf5d17253290a83 Mon Sep 17 00:00:00 2001 From: Mark Olesen Date: Mon, 20 Dec 2010 10:37:17 +0100 Subject: [PATCH] COMP: avoid ambiguous construct from tmp - utils/ pre+post processing --- .../miscellaneous/postChannel/collapse.H | 44 ++++++++++--------- .../createTurbulenceFields.C | 6 +-- .../velocityField/Lambda2/Lambda2.C | 6 ++- .../postProcessing/velocityField/Mach/Mach.C | 4 +- .../postProcessing/velocityField/Q/Q.C | 10 ++--- .../velocityField/flowType/flowType.C | 6 +-- .../wall/wallHeatFlux/wallHeatFlux.C | 6 ++- .../wall/yPlusLES/createFields.H | 2 +- .../postProcessing/wall/yPlusLES/yPlusLES.C | 2 +- .../applyBoundaryLayer/applyBoundaryLayer.C | 2 +- .../applyBoundaryLayer/createFields.H | 2 +- 11 files changed, 48 insertions(+), 42 deletions(-) diff --git a/applications/utilities/postProcessing/miscellaneous/postChannel/collapse.H b/applications/utilities/postProcessing/miscellaneous/postChannel/collapse.H index b3bf559411..4b042ab63f 100644 --- a/applications/utilities/postProcessing/miscellaneous/postChannel/collapse.H +++ b/applications/utilities/postProcessing/miscellaneous/postChannel/collapse.H @@ -1,41 +1,43 @@ - scalarField UMeanXvalues = channelIndexing.collapse + scalarField UMeanXvalues ( - UMean.component(vector::X)() + channelIndexing.collapse(UMean.component(vector::X)()) ); - scalarField UMeanYvalues = channelIndexing.collapse + scalarField UMeanYvalues ( - UMean.component(vector::Y)() + channelIndexing.collapse(UMean.component(vector::Y)()) ); - scalarField UMeanZvalues = channelIndexing.collapse + scalarField UMeanZvalues ( - UMean.component(vector::Z)() + channelIndexing.collapse(UMean.component(vector::Z)()) ); - scalarField RxxValues = channelIndexing.collapse(Rxx); - scalarField RyyValues = channelIndexing.collapse(Ryy); - scalarField RzzValues = channelIndexing.collapse(Rzz); - scalarField RxyValues = channelIndexing.collapse(Rxy, true); + scalarField RxxValues(channelIndexing.collapse(Rxx)); + scalarField RyyValues(channelIndexing.collapse(Ryy)); + scalarField RzzValues(channelIndexing.collapse(Rzz)); + scalarField RxyValues(channelIndexing.collapse(Rxy, true)); - scalarField pPrime2MeanValues = channelIndexing.collapse(pPrime2Mean); + scalarField pPrime2MeanValues(channelIndexing.collapse(pPrime2Mean)); /* - scalarField epsilonValues = channelIndexing.collapse(epsilonMean); + scalarField epsilonValues(channelIndexing.collapse(epsilonMean)); - scalarField nuMeanValues = channelIndexing.collapse(nuMean); - scalarField nuPrimeValues = channelIndexing.collapse(nuPrime); + scalarField nuMeanValues(channelIndexing.collapse(nuMean)); + scalarField nuPrimeValues(channelIndexing.collapse(nuPrime)); - scalarField gammaDotMeanValues = channelIndexing.collapse(gammaDotMean); - scalarField gammaDotPrimeValues = channelIndexing.collapse(gammaDotPrime); + scalarField gammaDotMeanValues(channelIndexing.collapse(gammaDotMean)); + scalarField gammaDotPrimeValues(channelIndexing.collapse(gammaDotPrime)); */ - scalarField urmsValues = sqrt(mag(RxxValues)); - scalarField vrmsValues = sqrt(mag(RyyValues)); - scalarField wrmsValues = sqrt(mag(RzzValues)); + scalarField urmsValues(sqrt(mag(RxxValues))); + scalarField vrmsValues(sqrt(mag(RyyValues))); + scalarField wrmsValues(sqrt(mag(RzzValues))); - scalarField kValues = - 0.5*(sqr(urmsValues) + sqr(vrmsValues) + sqr(wrmsValues)); + scalarField kValues + ( + 0.5*(sqr(urmsValues) + sqr(vrmsValues) + sqr(wrmsValues)) + ); const scalarField& y = channelIndexing.y(); diff --git a/applications/utilities/postProcessing/turbulence/createTurbulenceFields/createTurbulenceFields.C b/applications/utilities/postProcessing/turbulence/createTurbulenceFields/createTurbulenceFields.C index 761ea3bc65..fdd66e10e3 100644 --- a/applications/utilities/postProcessing/turbulence/createTurbulenceFields/createTurbulenceFields.C +++ b/applications/utilities/postProcessing/turbulence/createTurbulenceFields/createTurbulenceFields.C @@ -61,13 +61,13 @@ int main(int argc, char *argv[]) // Cache the turbulence fields Info<< "\nRetrieving field k from turbulence model" << endl; - const volScalarField k = RASModel->k(); + const volScalarField k(RASModel->k()); Info<< "\nRetrieving field epsilon from turbulence model" << endl; - const volScalarField epsilon = RASModel->epsilon(); + const volScalarField epsilon(RASModel->epsilon()); Info<< "\nRetrieving field R from turbulence model" << endl; - const volSymmTensorField R = RASModel->R(); + const volSymmTensorField R(RASModel->R()); // Check availability of tubulence fields diff --git a/applications/utilities/postProcessing/velocityField/Lambda2/Lambda2.C b/applications/utilities/postProcessing/velocityField/Lambda2/Lambda2.C index c089d48539..00e72e058d 100644 --- a/applications/utilities/postProcessing/velocityField/Lambda2/Lambda2.C +++ b/applications/utilities/postProcessing/velocityField/Lambda2/Lambda2.C @@ -55,8 +55,10 @@ void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh) const volTensorField gradU(fvc::grad(U)); - volTensorField SSplusWW = - (symm(gradU) & symm(gradU)) + (skew(gradU) & skew(gradU)); + volTensorField SSplusWW + ( + (symm(gradU) & symm(gradU)) + (skew(gradU) & skew(gradU)) + ); volScalarField Lambda2 ( diff --git a/applications/utilities/postProcessing/velocityField/Mach/Mach.C b/applications/utilities/postProcessing/velocityField/Mach/Mach.C index 304f30cdee..66f3d2c6f5 100644 --- a/applications/utilities/postProcessing/velocityField/Mach/Mach.C +++ b/applications/utilities/postProcessing/velocityField/Mach/Mach.C @@ -71,8 +71,8 @@ void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh) basicPsiThermo::New(mesh) ); - volScalarField Cp = thermo->Cp(); - volScalarField Cv = thermo->Cv(); + volScalarField Cp(thermo->Cp()); + volScalarField Cv(thermo->Cv()); MachPtr.set ( diff --git a/applications/utilities/postProcessing/velocityField/Q/Q.C b/applications/utilities/postProcessing/velocityField/Q/Q.C index 5e148a8041..bbd3f3422e 100644 --- a/applications/utilities/postProcessing/velocityField/Q/Q.C +++ b/applications/utilities/postProcessing/velocityField/Q/Q.C @@ -53,7 +53,7 @@ void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh) { Info<< " Reading U" << endl; volVectorField U(Uheader, mesh); - volTensorField gradU = fvc::grad(U); + volTensorField gradU(fvc::grad(U)); volScalarField Q ( @@ -72,11 +72,11 @@ void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh) // This is a second way of calculating Q, that delivers results // very close, but not identical to the first approach. - volSymmTensorField S = symm(gradU); // symmetric part of tensor - volTensorField W = skew(gradU); // anti-symmetric part + volSymmTensorField S(symm(gradU)); // symmetric part of tensor + volTensorField W(skew(gradU)); // anti-symmetric part - volScalarField SS = S&&S; - volScalarField WW = W&&W; + volScalarField SS(S && S); + volScalarField WW(W && W); volScalarField Q ( diff --git a/applications/utilities/postProcessing/velocityField/flowType/flowType.C b/applications/utilities/postProcessing/velocityField/flowType/flowType.C index 2a5936f0a0..bbee6e2be5 100644 --- a/applications/utilities/postProcessing/velocityField/flowType/flowType.C +++ b/applications/utilities/postProcessing/velocityField/flowType/flowType.C @@ -62,9 +62,9 @@ void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh) Info<< " Reading U" << endl; volVectorField U(Uheader, mesh); - volTensorField gradU = fvc::grad(U); - volScalarField magD = mag(symm(gradU)); - volScalarField magOmega = mag(skew(gradU)); + volTensorField gradU(fvc::grad(U)); + volScalarField magD(mag(symm(gradU))); + volScalarField magOmega (mag(skew(gradU))); dimensionedScalar smallMagD("smallMagD", magD.dimensions(), SMALL); Info<< " Calculating flowType" << endl; diff --git a/applications/utilities/postProcessing/wall/wallHeatFlux/wallHeatFlux.C b/applications/utilities/postProcessing/wall/wallHeatFlux/wallHeatFlux.C index c1bd9f05f4..f8f4a7c351 100644 --- a/applications/utilities/postProcessing/wall/wallHeatFlux/wallHeatFlux.C +++ b/applications/utilities/postProcessing/wall/wallHeatFlux/wallHeatFlux.C @@ -54,8 +54,10 @@ int main(int argc, char *argv[]) #include "createFields.H" - surfaceScalarField heatFlux = - fvc::interpolate(RASModel->alphaEff())*fvc::snGrad(h); + surfaceScalarField heatFlux + ( + fvc::interpolate(RASModel->alphaEff())*fvc::snGrad(h) + ); const surfaceScalarField::GeometricBoundaryField& patchHeatFlux = heatFlux.boundaryField(); diff --git a/applications/utilities/postProcessing/wall/yPlusLES/createFields.H b/applications/utilities/postProcessing/wall/yPlusLES/createFields.H index 83fb26a1d7..108aa9689f 100644 --- a/applications/utilities/postProcessing/wall/yPlusLES/createFields.H +++ b/applications/utilities/postProcessing/wall/yPlusLES/createFields.H @@ -21,4 +21,4 @@ autoPtr sgsModel incompressible::LESModel::New(U, phi, laminarTransport) ); -volScalarField::GeometricBoundaryField d = nearWallDist(mesh).y(); +volScalarField::GeometricBoundaryField d(nearWallDist(mesh).y()); diff --git a/applications/utilities/postProcessing/wall/yPlusLES/yPlusLES.C b/applications/utilities/postProcessing/wall/yPlusLES/yPlusLES.C index b104a114d5..2259aa0154 100644 --- a/applications/utilities/postProcessing/wall/yPlusLES/yPlusLES.C +++ b/applications/utilities/postProcessing/wall/yPlusLES/yPlusLES.C @@ -100,7 +100,7 @@ int main(int argc, char *argv[]) ); volScalarField::GeometricBoundaryField d = nearWallDist(mesh).y(); - volScalarField nuEff = sgsModel->nuEff(); + volScalarField nuEff(sgsModel->nuEff()); const fvPatchList& patches = mesh.boundary(); diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C index 8fe46ba299..5b105742fa 100644 --- a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C +++ b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C @@ -104,7 +104,7 @@ int main(int argc, char *argv[]) // Calculate nut tmp tnut = turbulence->nut(); volScalarField& nut = tnut(); - volScalarField S = mag(dev(symm(fvc::grad(U)))); + volScalarField S(mag(dev(symm(fvc::grad(U))))); nut = sqr(kappa*min(y, ybl))*::sqrt(2)*S; if (args.optionFound("writenut")) diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/createFields.H b/applications/utilities/preProcessing/applyBoundaryLayer/createFields.H index b8f5e0dade..7447ad38e8 100644 --- a/applications/utilities/preProcessing/applyBoundaryLayer/createFields.H +++ b/applications/utilities/preProcessing/applyBoundaryLayer/createFields.H @@ -47,7 +47,7 @@ License ); Info<< "Calculating wall distance field" << endl; - volScalarField y = wallDist(mesh).y(); + volScalarField y(wallDist(mesh).y()); // Set the mean boundary-layer thickness dimensionedScalar ybl("ybl", dimLength, 0);