diff --git a/applications/solvers/multiphase/bubbleFoam/UEqns.H b/applications/solvers/multiphase/bubbleFoam/UEqns.H index 64cc5db462..82d53d84b0 100644 --- a/applications/solvers/multiphase/bubbleFoam/UEqns.H +++ b/applications/solvers/multiphase/bubbleFoam/UEqns.H @@ -2,12 +2,14 @@ fvVectorMatrix UaEqn(Ua, Ua.dimensions()*dimVol/dimTime); fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime); { - volTensorField Rca = -nuEffa*(fvc::grad(Ua)().T()); + volTensorField Rca(-nuEffa*(fvc::grad(Ua)().T())); Rca = Rca + (2.0/3.0)*sqr(Ct)*I*k - (2.0/3.0)*I*tr(Rca); - surfaceScalarField phiRa = + surfaceScalarField phiRa + ( - fvc::interpolate(nuEffa) - *mesh.magSf()*fvc::snGrad(alpha)/fvc::interpolate(alpha + scalar(0.001)); + *mesh.magSf()*fvc::snGrad(alpha)/fvc::interpolate(alpha + scalar(0.001)) + ); UaEqn = ( @@ -34,12 +36,14 @@ fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime); UaEqn.relax(); - volTensorField Rcb = -nuEffb*fvc::grad(Ub)().T(); + volTensorField Rcb(-nuEffb*fvc::grad(Ub)().T()); Rcb = Rcb + (2.0/3.0)*I*k - (2.0/3.0)*I*tr(Rcb); - surfaceScalarField phiRb = + surfaceScalarField phiRb + ( - fvc::interpolate(nuEffb) - *mesh.magSf()*fvc::snGrad(beta)/fvc::interpolate(beta + scalar(0.001)); + *mesh.magSf()*fvc::snGrad(beta)/fvc::interpolate(beta + scalar(0.001)) + ); UbEqn = ( diff --git a/applications/solvers/multiphase/bubbleFoam/alphaEqn.H b/applications/solvers/multiphase/bubbleFoam/alphaEqn.H index db3b595625..772008bc89 100644 --- a/applications/solvers/multiphase/bubbleFoam/alphaEqn.H +++ b/applications/solvers/multiphase/bubbleFoam/alphaEqn.H @@ -1,7 +1,7 @@ { word scheme("div(phi,alpha)"); - surfaceScalarField phir = phia - phib; + surfaceScalarField phir(phia - phib); Info<< "Max Ur Courant Number = " << ( diff --git a/applications/solvers/multiphase/bubbleFoam/createFields.H b/applications/solvers/multiphase/bubbleFoam/createFields.H index 57bb7e88d5..efb41613a6 100644 --- a/applications/solvers/multiphase/bubbleFoam/createFields.H +++ b/applications/solvers/multiphase/bubbleFoam/createFields.H @@ -171,15 +171,19 @@ Info<< "Calculating field DDtUa and DDtUb\n" << endl; - volVectorField DDtUa = + volVectorField DDtUa + ( fvc::ddt(Ua) + fvc::div(phia, Ua) - - fvc::div(phia)*Ua; + - fvc::div(phia)*Ua + ); - volVectorField DDtUb = + volVectorField DDtUb + ( fvc::ddt(Ub) + fvc::div(phib, Ub) - - fvc::div(phib)*Ub; + - fvc::div(phib)*Ub + ); Info<< "Calculating field g.h\n" << endl; diff --git a/applications/solvers/multiphase/bubbleFoam/kEpsilon.H b/applications/solvers/multiphase/bubbleFoam/kEpsilon.H index 84dadd7acb..309d4eba80 100644 --- a/applications/solvers/multiphase/bubbleFoam/kEpsilon.H +++ b/applications/solvers/multiphase/bubbleFoam/kEpsilon.H @@ -6,7 +6,7 @@ if (turbulence) } tmp tgradUb = fvc::grad(Ub); - volScalarField G = 2*nutb*(tgradUb() && dev(symm(tgradUb()))); + volScalarField G(2*nutb*(tgradUb() && dev(symm(tgradUb())))); tgradUb.clear(); #include "wallFunctions.H" diff --git a/applications/solvers/multiphase/bubbleFoam/liftDragCoeffs.H b/applications/solvers/multiphase/bubbleFoam/liftDragCoeffs.H index f9ae6ab2ca..5d00e39820 100644 --- a/applications/solvers/multiphase/bubbleFoam/liftDragCoeffs.H +++ b/applications/solvers/multiphase/bubbleFoam/liftDragCoeffs.H @@ -1,11 +1,15 @@ - volVectorField Ur = Ua - Ub; - volScalarField magUr = mag(Ur); + volVectorField Ur(Ua - Ub); + volScalarField magUr(mag(Ur)); - volScalarField CdaMagUr = - (24.0*nub/da)*(scalar(1) + 0.15*pow(da*magUr/nub, 0.687)); + volScalarField CdaMagUr + ( + (24.0*nub/da)*(scalar(1) + 0.15*pow(da*magUr/nub, 0.687)) + ); - volScalarField CdbMagUr = - (24.0*nua/db)*(scalar(1) + 0.15*pow(db*magUr/nua, 0.687)); + volScalarField CdbMagUr + ( + (24.0*nua/db)*(scalar(1) + 0.15*pow(db*magUr/nua, 0.687)) + ); volScalarField dragCoef ( @@ -13,4 +17,7 @@ 0.75*(beta*rhob*CdaMagUr/da + alpha*rhoa*CdbMagUr/db) ); - volVectorField liftCoeff = Cl*(beta*rhob + alpha*rhoa)*(Ur ^ fvc::curl(U)); + volVectorField liftCoeff + ( + Cl*(beta*rhob + alpha*rhoa)*(Ur ^ fvc::curl(U)) + ); diff --git a/applications/solvers/multiphase/bubbleFoam/pEqn.H b/applications/solvers/multiphase/bubbleFoam/pEqn.H index 8e257c8fa9..faae92d397 100644 --- a/applications/solvers/multiphase/bubbleFoam/pEqn.H +++ b/applications/solvers/multiphase/bubbleFoam/pEqn.H @@ -1,20 +1,24 @@ { - surfaceScalarField alphaf = fvc::interpolate(alpha); - surfaceScalarField betaf = scalar(1) - alphaf; + surfaceScalarField alphaf(fvc::interpolate(alpha)); + surfaceScalarField betaf(scalar(1) - alphaf); - volScalarField rUaA = 1.0/UaEqn.A(); - volScalarField rUbA = 1.0/UbEqn.A(); + volScalarField rUaA(1.0/UaEqn.A()); + volScalarField rUbA(1.0/UbEqn.A()); - surfaceScalarField rUaAf = fvc::interpolate(rUaA); - surfaceScalarField rUbAf = fvc::interpolate(rUbA); + surfaceScalarField rUaAf(fvc::interpolate(rUaA)); + surfaceScalarField rUbAf(fvc::interpolate(rUbA)); Ua = rUaA*UaEqn.H(); Ub = rUbA*UbEqn.H(); - surfaceScalarField phiDraga = - fvc::interpolate(beta/rhoa*dragCoef*rUaA)*phib + rUaAf*(g & mesh.Sf()); - surfaceScalarField phiDragb = - fvc::interpolate(alpha/rhob*dragCoef*rUbA)*phia + rUbAf*(g & mesh.Sf()); + surfaceScalarField phiDraga + ( + fvc::interpolate(beta/rhoa*dragCoef*rUaA)*phib + rUaAf*(g & mesh.Sf()) + ); + surfaceScalarField phiDragb + ( + fvc::interpolate(alpha/rhob*dragCoef*rUbA)*phia + rUbAf*(g & mesh.Sf()) + ); forAll(p.boundaryField(), patchi) { @@ -50,7 +54,7 @@ if (nonOrth == nNonOrthCorr) { - surfaceScalarField SfGradp = pEqn.flux()/Dp; + surfaceScalarField SfGradp(pEqn.flux()/Dp); phia -= rUaAf*SfGradp/rhoa; phib -= rUbAf*SfGradp/rhob; diff --git a/applications/solvers/multiphase/bubbleFoam/wallFunctions.H b/applications/solvers/multiphase/bubbleFoam/wallFunctions.H index fe5d16909e..a998d4c3dd 100644 --- a/applications/solvers/multiphase/bubbleFoam/wallFunctions.H +++ b/applications/solvers/multiphase/bubbleFoam/wallFunctions.H @@ -34,7 +34,7 @@ { const scalarField& nuw = nutb.boundaryField()[patchi]; - scalarField magFaceGradU = mag(U.boundaryField()[patchi].snGrad()); + scalarField magFaceGradU(mag(U.boundaryField()[patchi].snGrad())); forAll(currPatch, facei) { diff --git a/applications/solvers/multiphase/cavitatingFoam/CourantNo.H b/applications/solvers/multiphase/cavitatingFoam/CourantNo.H index 8b26ba033e..cb72ce18f7 100644 --- a/applications/solvers/multiphase/cavitatingFoam/CourantNo.H +++ b/applications/solvers/multiphase/cavitatingFoam/CourantNo.H @@ -35,8 +35,10 @@ scalar acousticCoNum = 0.0; if (mesh.nInternalFaces()) { - scalarField sumPhi = - fvc::surfaceSum(mag(phiv))().internalField(); + scalarField sumPhi + ( + fvc::surfaceSum(mag(phiv))().internalField() + ); CoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue(); diff --git a/applications/solvers/multiphase/cavitatingFoam/pEqn.H b/applications/solvers/multiphase/cavitatingFoam/pEqn.H index e6606b7b62..6b6fc73eb9 100644 --- a/applications/solvers/multiphase/cavitatingFoam/pEqn.H +++ b/applications/solvers/multiphase/cavitatingFoam/pEqn.H @@ -9,18 +9,18 @@ )/psi; } - surfaceScalarField rhof = fvc::interpolate(rho, "rhof"); + surfaceScalarField rhof(fvc::interpolate(rho, "rhof")); - volScalarField rAU = 1.0/UEqn.A(); + volScalarField rAU(1.0/UEqn.A()); surfaceScalarField rAUf("rAUf", rhof*fvc::interpolate(rAU)); - volVectorField HbyA = rAU*UEqn.H(); + volVectorField HbyA(rAU*UEqn.H()); phiv = (fvc::interpolate(HbyA) & mesh.Sf()) + fvc::ddtPhiCorr(rAU, rho, U, phiv); p.boundaryField().updateCoeffs(); - surfaceScalarField phiGradp = rAUf*mesh.magSf()*fvc::snGrad(p); + surfaceScalarField phiGradp(rAUf*mesh.magSf()*fvc::snGrad(p)); phiv -= phiGradp/rhof; diff --git a/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H b/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H index f545262767..f5717a654d 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H +++ b/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H @@ -2,7 +2,7 @@ word alphaScheme("div(phi,alpha)"); word alpharScheme("div(phirb,alpha)"); - surfaceScalarField phir = phic*interface.nHatf(); + surfaceScalarField phir(phic*interface.nHatf()); for (int gCorr=0; gCorr 1) { dimensionedScalar totalDeltaT = runTime.deltaT(); - surfaceScalarField rhoPhiSum = 0.0*rhoPhi; + surfaceScalarField rhoPhiSum(0.0*rhoPhi); for ( diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/alphaEqnsSubCycle.H b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/alphaEqnsSubCycle.H index e668296b6d..a973c23dfe 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/alphaEqnsSubCycle.H +++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/alphaEqnsSubCycle.H @@ -9,17 +9,17 @@ readLabel(piso.lookup("nAlphaSubCycles")) ); - surfaceScalarField phic = mag(phi/mesh.magSf()); + surfaceScalarField phic(mag(phi/mesh.magSf())); phic = min(interface.cAlpha()*phic, max(phic)); fvc::makeAbsolute(phi, U); - volScalarField divU = fvc::div(phi); + volScalarField divU(fvc::div(phi)); fvc::makeRelative(phi, U); if (nAlphaSubCycles > 1) { dimensionedScalar totalDeltaT = runTime.deltaT(); - surfaceScalarField rhoPhiSum = 0.0*rhoPhi; + surfaceScalarField rhoPhiSum(0.0*rhoPhi); for ( diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C index c2f1d04706..8519d23f67 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C +++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C @@ -79,7 +79,7 @@ int main(int argc, char *argv[]) { // Store divU from the previous mesh for the correctPhi - volScalarField divU = fvc::div(phi); + volScalarField divU(fvc::div(phi)); scalar timeBeforeMeshUpdate = runTime.elapsedCpuTime(); diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/pEqn.H b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/pEqn.H index ace51b5d84..0c1ad0192d 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/pEqn.H +++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); tmp p_rghEqnComp; diff --git a/applications/solvers/multiphase/compressibleInterFoam/createFields.H b/applications/solvers/multiphase/compressibleInterFoam/createFields.H index fab0df9284..c598cb75ce 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/createFields.H +++ b/applications/solvers/multiphase/compressibleInterFoam/createFields.H @@ -105,8 +105,8 @@ ) ); - volScalarField rho1 = rho10 + psi1*p; - volScalarField rho2 = rho20 + psi2*p; + volScalarField rho1(rho10 + psi1*p); + volScalarField rho2(rho20 + psi2*p); volScalarField rho ( @@ -138,8 +138,10 @@ fvc::interpolate(rho)*phi ); - volScalarField dgdt = - pos(alpha2)*fvc::div(phi)/max(alpha2, scalar(0.0001)); + volScalarField dgdt + ( + pos(alpha2)*fvc::div(phi)/max(alpha2, scalar(0.0001)) + ); // Construct interface from alpha1 distribution interfaceProperties interface(alpha1, U, twoPhaseProperties); diff --git a/applications/solvers/multiphase/compressibleInterFoam/pEqn.H b/applications/solvers/multiphase/compressibleInterFoam/pEqn.H index 7054fe1d44..047f2f3c4c 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/pEqn.H +++ b/applications/solvers/multiphase/compressibleInterFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); tmp p_rghEqnComp; diff --git a/applications/solvers/multiphase/interFoam/LTSInterFoam/MULESTemplates.C b/applications/solvers/multiphase/interFoam/LTSInterFoam/MULESTemplates.C index 280d7d6d55..d6ab8bfc3e 100644 --- a/applications/solvers/multiphase/interFoam/LTSInterFoam/MULESTemplates.C +++ b/applications/solvers/multiphase/interFoam/LTSInterFoam/MULESTemplates.C @@ -56,7 +56,7 @@ void Foam::MULES::explicitLTSSolve const fvMesh& mesh = psi.mesh(); psi.correctBoundaryConditions(); - surfaceScalarField phiBD = upwind(psi.mesh(), phi).flux(psi); + surfaceScalarField phiBD(upwind(psi.mesh(), phi).flux(psi)); surfaceScalarField& phiCorr = phiPsi; phiCorr -= phiBD; @@ -168,9 +168,11 @@ void Foam::MULES::implicitSolve scalarField allCoLambda(mesh.nFaces()); { - surfaceScalarField Cof = + surfaceScalarField Cof + ( mesh.time().deltaT()*mesh.surfaceInterpolation::deltaCoeffs() - *mag(phi)/mesh.magSf(); + *mag(phi)/mesh.magSf() + ); slicedSurfaceScalarField CoLambda ( @@ -226,7 +228,7 @@ void Foam::MULES::implicitSolve - Su ); - surfaceScalarField phiBD = psiConvectionDiffusion.flux(); + surfaceScalarField phiBD(psiConvectionDiffusion.flux()); surfaceScalarField& phiCorr = phiPsi; phiCorr -= phiBD; @@ -408,7 +410,7 @@ void Foam::MULES::limiter if (psiPf.coupled()) { - scalarField psiPNf = psiPf.patchNeighbourField(); + scalarField psiPNf(psiPf.patchNeighbourField()); forAll(phiCorrPf, pFacei) { diff --git a/applications/solvers/multiphase/interFoam/LTSInterFoam/alphaEqn.H b/applications/solvers/multiphase/interFoam/LTSInterFoam/alphaEqn.H index a4cb6e7688..0c2cf71e4e 100644 --- a/applications/solvers/multiphase/interFoam/LTSInterFoam/alphaEqn.H +++ b/applications/solvers/multiphase/interFoam/LTSInterFoam/alphaEqn.H @@ -2,13 +2,14 @@ word alphaScheme("div(phi,alpha)"); word alpharScheme("div(phirb,alpha)"); - surfaceScalarField phic = mag(phi/mesh.magSf()); + surfaceScalarField phic(mag(phi/mesh.magSf())); phic = min(interface.cAlpha()*phic, max(phic)); - surfaceScalarField phir = phic*interface.nHatf(); + surfaceScalarField phir(phic*interface.nHatf()); for (int aCorr=0; aCorr 1) { dimensionedScalar totalDeltaT = runTime.deltaT(); - surfaceScalarField rhoPhiSum = 0.0*rhoPhi; + surfaceScalarField rhoPhiSum(0.0*rhoPhi); for ( diff --git a/applications/solvers/multiphase/interFoam/MRFInterFoam/pEqn.H b/applications/solvers/multiphase/interFoam/MRFInterFoam/pEqn.H index 7944df4814..e4f8247bac 100644 --- a/applications/solvers/multiphase/interFoam/MRFInterFoam/pEqn.H +++ b/applications/solvers/multiphase/interFoam/MRFInterFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); U = rAU*UEqn.H(); surfaceScalarField phiU diff --git a/applications/solvers/multiphase/interFoam/alphaCourantNo.H b/applications/solvers/multiphase/interFoam/alphaCourantNo.H index e5edfc76ca..d122753d46 100644 --- a/applications/solvers/multiphase/interFoam/alphaCourantNo.H +++ b/applications/solvers/multiphase/interFoam/alphaCourantNo.H @@ -39,9 +39,11 @@ scalar meanAlphaCoNum = 0.0; if (mesh.nInternalFaces()) { - scalarField sumPhi = + scalarField sumPhi + ( pos(alpha1 - 0.01)*pos(0.99 - alpha1) - *fvc::surfaceSum(mag(phi))().internalField(); + *fvc::surfaceSum(mag(phi))().internalField() + ); alphaCoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue(); diff --git a/applications/solvers/multiphase/interFoam/alphaEqn.H b/applications/solvers/multiphase/interFoam/alphaEqn.H index 0b2fb4ebf8..d19ae9bb70 100644 --- a/applications/solvers/multiphase/interFoam/alphaEqn.H +++ b/applications/solvers/multiphase/interFoam/alphaEqn.H @@ -2,13 +2,14 @@ word alphaScheme("div(phi,alpha)"); word alpharScheme("div(phirb,alpha)"); - surfaceScalarField phic = mag(phi/mesh.magSf()); + surfaceScalarField phic(mag(phi/mesh.magSf())); phic = min(interface.cAlpha()*phic, max(phic)); - surfaceScalarField phir = phic*interface.nHatf(); + surfaceScalarField phir(phic*interface.nHatf()); for (int aCorr=0; aCorr 1) { dimensionedScalar totalDeltaT = runTime.deltaT(); - surfaceScalarField rhoPhiSum = 0.0*rhoPhi; + surfaceScalarField rhoPhiSum(0.0*rhoPhi); for ( diff --git a/applications/solvers/multiphase/interFoam/interDyMFoam/pEqn.H b/applications/solvers/multiphase/interFoam/interDyMFoam/pEqn.H index 59c0cd1929..5bf3ca9418 100644 --- a/applications/solvers/multiphase/interFoam/interDyMFoam/pEqn.H +++ b/applications/solvers/multiphase/interFoam/interDyMFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); U = rAU*UEqn.H(); surfaceScalarField phiU("phiU", (fvc::interpolate(U) & mesh.Sf())); diff --git a/applications/solvers/multiphase/interFoam/pEqn.H b/applications/solvers/multiphase/interFoam/pEqn.H index 6437572c3b..db1eccc5a5 100644 --- a/applications/solvers/multiphase/interFoam/pEqn.H +++ b/applications/solvers/multiphase/interFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); U = rAU*UEqn.H(); surfaceScalarField phiU diff --git a/applications/solvers/multiphase/interMixingFoam/alphaCourantNo.H b/applications/solvers/multiphase/interMixingFoam/alphaCourantNo.H index 105149f87e..821064fae1 100644 --- a/applications/solvers/multiphase/interMixingFoam/alphaCourantNo.H +++ b/applications/solvers/multiphase/interMixingFoam/alphaCourantNo.H @@ -39,11 +39,14 @@ scalar meanAlphaCoNum = 0.0; if (mesh.nInternalFaces()) { - scalarField sumPhi = max + scalarField sumPhi ( - pos(alpha1 - 0.01)*pos(0.99 - alpha1), - pos(alpha2 - 0.01)*pos(0.99 - alpha2) - )*fvc::surfaceSum(mag(phi))().internalField(); + max + ( + pos(alpha1 - 0.01)*pos(0.99 - alpha1), + pos(alpha2 - 0.01)*pos(0.99 - alpha2) + )*fvc::surfaceSum(mag(phi))().internalField() + ); alphaCoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue(); diff --git a/applications/solvers/multiphase/interMixingFoam/alphaEqns.H b/applications/solvers/multiphase/interMixingFoam/alphaEqns.H index 00295354d4..15759adfe4 100644 --- a/applications/solvers/multiphase/interMixingFoam/alphaEqns.H +++ b/applications/solvers/multiphase/interMixingFoam/alphaEqns.H @@ -40,7 +40,8 @@ // Create the complete convection flux for alpha1 - surfaceScalarField phiAlpha1 = + surfaceScalarField phiAlpha1 + ( fvc::flux ( phi, @@ -58,11 +59,14 @@ -fvc::flux(-phir, alpha3, alpharScheme), alpha1, alpharScheme - ); + ) + ); // Create the bounded (upwind) flux for alpha1 - surfaceScalarField phiAlpha1BD = - upwind(mesh, phi).flux(alpha1); + surfaceScalarField phiAlpha1BD + ( + upwind(mesh, phi).flux(alpha1) + ); // Calculate the flux correction for alpha1 phiAlpha1 -= phiAlpha1BD; @@ -83,7 +87,8 @@ ); // Create the complete flux for alpha2 - surfaceScalarField phiAlpha2 = + surfaceScalarField phiAlpha2 + ( fvc::flux ( phi, @@ -95,11 +100,14 @@ -fvc::flux(phir, alpha1, alpharScheme), alpha2, alpharScheme - ); + ) + ); // Create the bounded (upwind) flux for alpha2 - surfaceScalarField phiAlpha2BD = - upwind(mesh, phi).flux(alpha2); + surfaceScalarField phiAlpha2BD + ( + upwind(mesh, phi).flux(alpha2) + ); // Calculate the flux correction for alpha2 phiAlpha2 -= phiAlpha2BD; @@ -127,8 +135,8 @@ solve(fvm::ddt(alpha1) + fvc::div(phiAlpha1)); // Create the diffusion coefficients for alpha2<->alpha3 - volScalarField Dc23 = D23*max(alpha3, scalar(0))*pos(alpha2); - volScalarField Dc32 = D23*max(alpha2, scalar(0))*pos(alpha3); + volScalarField Dc23(D23*max(alpha3, scalar(0))*pos(alpha2)); + volScalarField Dc32(D23*max(alpha2, scalar(0))*pos(alpha3)); // Add the diffusive flux for alpha3->alpha2 phiAlpha2 -= fvc::interpolate(Dc32)*mesh.magSf()*fvc::snGrad(alpha1); diff --git a/applications/solvers/multiphase/interMixingFoam/alphaEqnsSubCycle.H b/applications/solvers/multiphase/interMixingFoam/alphaEqnsSubCycle.H index 765087a183..81e09f9b90 100644 --- a/applications/solvers/multiphase/interMixingFoam/alphaEqnsSubCycle.H +++ b/applications/solvers/multiphase/interMixingFoam/alphaEqnsSubCycle.H @@ -10,7 +10,7 @@ label nAlphaSubCycles if (nAlphaSubCycles > 1) { - surfaceScalarField rhoPhiSum = 0.0*rhoPhi; + surfaceScalarField rhoPhiSum(0.0*rhoPhi); dimensionedScalar totalDeltaT = runTime.deltaT(); for @@ -33,7 +33,7 @@ else interface.correct(); { - volScalarField rhoNew = alpha1*rho1 + alpha2*rho2 + alpha3*rho3; + volScalarField rhoNew(alpha1*rho1 + alpha2*rho2 + alpha3*rho3); //solve(fvm::ddt(rho) + fvc::div(rhoPhi)); //Info<< "density error = " diff --git a/applications/solvers/multiphase/interMixingFoam/incompressibleThreePhaseMixture/threePhaseMixture.C b/applications/solvers/multiphase/interMixingFoam/incompressibleThreePhaseMixture/threePhaseMixture.C index 4ea4c4d989..163d5a5410 100644 --- a/applications/solvers/multiphase/interMixingFoam/incompressibleThreePhaseMixture/threePhaseMixture.C +++ b/applications/solvers/multiphase/interMixingFoam/incompressibleThreePhaseMixture/threePhaseMixture.C @@ -133,9 +133,9 @@ Foam::tmp Foam::threePhaseMixture::mu() const Foam::tmp Foam::threePhaseMixture::muf() const { - surfaceScalarField alpha1f = fvc::interpolate(alpha1_); - surfaceScalarField alpha2f = fvc::interpolate(alpha2_); - surfaceScalarField alpha3f = fvc::interpolate(alpha3_); + surfaceScalarField alpha1f(fvc::interpolate(alpha1_)); + surfaceScalarField alpha2f(fvc::interpolate(alpha2_)); + surfaceScalarField alpha3f(fvc::interpolate(alpha3_)); return tmp ( @@ -152,9 +152,9 @@ Foam::tmp Foam::threePhaseMixture::muf() const Foam::tmp Foam::threePhaseMixture::nuf() const { - surfaceScalarField alpha1f = fvc::interpolate(alpha1_); - surfaceScalarField alpha2f = fvc::interpolate(alpha2_); - surfaceScalarField alpha3f = fvc::interpolate(alpha3_); + surfaceScalarField alpha1f(fvc::interpolate(alpha1_)); + surfaceScalarField alpha2f(fvc::interpolate(alpha2_)); + surfaceScalarField alpha3f(fvc::interpolate(alpha3_)); return tmp ( diff --git a/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C b/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C index 48089e7c92..ad9e9fa1c5 100644 --- a/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C +++ b/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.C @@ -81,31 +81,35 @@ void Foam::threePhaseInterfaceProperties::correctContactAngle refCast (alpha3[patchi]); - scalarField twoPhaseAlpha2 = max(a2cap, scalar(0)); - scalarField twoPhaseAlpha3 = max(a3cap, scalar(0)); + scalarField twoPhaseAlpha2(max(a2cap, scalar(0))); + scalarField twoPhaseAlpha3(max(a3cap, scalar(0))); - scalarField sumTwoPhaseAlpha = - twoPhaseAlpha2 + twoPhaseAlpha3 + SMALL; + scalarField sumTwoPhaseAlpha + ( + twoPhaseAlpha2 + twoPhaseAlpha3 + SMALL + ); twoPhaseAlpha2 /= sumTwoPhaseAlpha; twoPhaseAlpha3 /= sumTwoPhaseAlpha; fvsPatchVectorField& nHatp = nHatb[patchi]; - scalarField theta = + scalarField theta + ( convertToRad - *( + * ( twoPhaseAlpha2*(180 - a2cap.theta(U[patchi], nHatp)) + twoPhaseAlpha3*(180 - a3cap.theta(U[patchi], nHatp)) - ); + ) + ); - vectorField nf = boundary[patchi].nf(); + vectorField nf(boundary[patchi].nf()); // Reset nHatPatch to correspond to the contact angle - scalarField a12 = nHatp & nf; + scalarField a12(nHatp & nf); - scalarField b1 = cos(theta); + scalarField b1(cos(theta)); scalarField b2(nHatp.size()); @@ -114,10 +118,10 @@ void Foam::threePhaseInterfaceProperties::correctContactAngle b2[facei] = cos(acos(a12[facei]) - theta[facei]); } - scalarField det = 1.0 - a12*a12; + scalarField det(1.0 - a12*a12); - scalarField a = (b1 - a12*b2)/det; - scalarField b = (b2 - a12*b1)/det; + scalarField a((b1 - a12*b2)/det); + scalarField b((b2 - a12*b1)/det); nHatp = a*nf + b*nHatp; @@ -135,13 +139,13 @@ void Foam::threePhaseInterfaceProperties::calculateK() const surfaceVectorField& Sf = mesh.Sf(); // Cell gradient of alpha - volVectorField gradAlpha = fvc::grad(alpha1); + volVectorField gradAlpha(fvc::grad(alpha1)); // Interpolated face-gradient of alpha - surfaceVectorField gradAlphaf = fvc::interpolate(gradAlpha); + surfaceVectorField gradAlphaf(fvc::interpolate(gradAlpha)); // Face unit interface normal - surfaceVectorField nHatfv = gradAlphaf/(mag(gradAlphaf) + deltaN_); + surfaceVectorField nHatfv(gradAlphaf/(mag(gradAlphaf) + deltaN_)); correctContactAngle(nHatfv.boundaryField()); // Face unit interface normal flux diff --git a/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.H b/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.H index d21409408f..d136aab5c0 100644 --- a/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.H +++ b/applications/solvers/multiphase/interMixingFoam/threePhaseInterfaceProperties/threePhaseInterfaceProperties.H @@ -125,8 +125,8 @@ public: tmp sigma() const { - volScalarField limitedAlpha2 = max(mixture_.alpha2(), scalar(0)); - volScalarField limitedAlpha3 = max(mixture_.alpha3(), scalar(0)); + volScalarField limitedAlpha2(max(mixture_.alpha2(), scalar(0))); + volScalarField limitedAlpha3(max(mixture_.alpha3(), scalar(0))); return (limitedAlpha2*sigma12_ + limitedAlpha3*sigma13_) diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H b/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H index 05226dc72f..9ce20db849 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H +++ b/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H @@ -6,7 +6,8 @@ for (int gCorr=0; gCorr > vDotAlphal = twoPhaseProperties->vDotAlphal(); diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqnSubCycle.H b/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqnSubCycle.H index df45bc6482..eea578e6e5 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqnSubCycle.H +++ b/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqnSubCycle.H @@ -21,10 +21,10 @@ surfaceScalarField rhoPhi readLabel(piso.lookup("nAlphaSubCycles")) ); - surfaceScalarField phic = mag(phi/mesh.magSf()); + surfaceScalarField phic(mag(phi/mesh.magSf())); phic = min(interface.cAlpha()*phic, max(phic)); - volScalarField divU = fvc::div(phi); + volScalarField divU(fvc::div(phi)); dimensionedScalar totalDeltaT = runTime.deltaT(); diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/pEqn.H b/applications/solvers/multiphase/interPhaseChangeFoam/pEqn.H index 4f4e68cb59..925e10aad2 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/pEqn.H +++ b/applications/solvers/multiphase/interPhaseChangeFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); U = rAU*UEqn.H(); diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.C b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.C index 4075122116..1459fdc2fc 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.C +++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.C @@ -68,7 +68,7 @@ Foam::Pair > Foam::phaseChangeTwoPhaseMixtures::Kunz::mDotAlphal() const { const volScalarField& p = alpha1_.db().lookupObject("p"); - volScalarField limitedAlpha1 = min(max(alpha1_, scalar(0)), scalar(1)); + volScalarField limitedAlpha1(min(max(alpha1_, scalar(0)), scalar(1))); return Pair > ( @@ -83,7 +83,7 @@ Foam::Pair > Foam::phaseChangeTwoPhaseMixtures::Kunz::mDotP() const { const volScalarField& p = alpha1_.db().lookupObject("p"); - volScalarField limitedAlpha1 = min(max(alpha1_, scalar(0)), scalar(1)); + volScalarField limitedAlpha1(min(max(alpha1_, scalar(0)), scalar(1))); return Pair > ( diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.C b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.C index 66c279df2c..250fb00578 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.C +++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.C @@ -80,7 +80,7 @@ Foam::Pair > Foam::phaseChangeTwoPhaseMixtures::Merkle::mDotP() const { const volScalarField& p = alpha1_.db().lookupObject("p"); - volScalarField limitedAlpha1 = min(max(alpha1_, scalar(0)), scalar(1)); + volScalarField limitedAlpha1(min(max(alpha1_, scalar(0)), scalar(1))); return Pair > ( diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.C b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.C index 61e57b3a5e..fb390277fb 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.C +++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.C @@ -97,9 +97,11 @@ Foam::phaseChangeTwoPhaseMixtures::SchnerrSauer::pCoeff const volScalarField& p ) const { - volScalarField limitedAlpha1 = min(max(alpha1_, scalar(0)), scalar(1)); - volScalarField rho = - (limitedAlpha1*rho1() + (scalar(1) - limitedAlpha1)*rho2()); + volScalarField limitedAlpha1(min(max(alpha1_, scalar(0)), scalar(1))); + volScalarField rho + ( + limitedAlpha1*rho1() + (scalar(1) - limitedAlpha1)*rho2() + ); return (3*rho1()*rho2())*sqrt(2/(3*rho1())) @@ -111,9 +113,9 @@ Foam::Pair > Foam::phaseChangeTwoPhaseMixtures::SchnerrSauer::mDotAlphal() const { const volScalarField& p = alpha1_.db().lookupObject("p"); - volScalarField limitedAlpha1 = min(max(alpha1_, scalar(0)), scalar(1)); + volScalarField limitedAlpha1(min(max(alpha1_, scalar(0)), scalar(1))); - volScalarField pCoeff = this->pCoeff(p); + volScalarField pCoeff(this->pCoeff(p)); return Pair > ( @@ -128,9 +130,9 @@ Foam::Pair > Foam::phaseChangeTwoPhaseMixtures::SchnerrSauer::mDotP() const { const volScalarField& p = alpha1_.db().lookupObject("p"); - volScalarField limitedAlpha1 = min(max(alpha1_, scalar(0)), scalar(1)); + volScalarField limitedAlpha1(min(max(alpha1_, scalar(0)), scalar(1))); - volScalarField apCoeff = limitedAlpha1*pCoeff(p); + volScalarField apCoeff(limitedAlpha1*pCoeff(p)); return Pair > ( diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.C b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.C index 0ae535dd51..374d61ddfa 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.C +++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.C @@ -54,7 +54,7 @@ Foam::phaseChangeTwoPhaseMixture::phaseChangeTwoPhaseMixture Foam::Pair > Foam::phaseChangeTwoPhaseMixture::vDotAlphal() const { - volScalarField alphalCoeff = 1.0/rho1() - alpha1_*(1.0/rho1() - 1.0/rho2()); + volScalarField alphalCoeff(1.0/rho1() - alpha1_*(1.0/rho1() - 1.0/rho2())); Pair > mDotAlphal = this->mDotAlphal(); return Pair > diff --git a/applications/solvers/multiphase/multiphaseInterFoam/MRFMultiphaseInterFoam/pEqn.H b/applications/solvers/multiphase/multiphaseInterFoam/MRFMultiphaseInterFoam/pEqn.H index 3f59b6b3fc..8be2e12bf6 100644 --- a/applications/solvers/multiphase/multiphaseInterFoam/MRFMultiphaseInterFoam/pEqn.H +++ b/applications/solvers/multiphase/multiphaseInterFoam/MRFMultiphaseInterFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); U = rAU*UEqn.H(); diff --git a/applications/solvers/multiphase/multiphaseInterFoam/alphaCourantNo.H b/applications/solvers/multiphase/multiphaseInterFoam/alphaCourantNo.H index a63ab62103..bfcb18f2c1 100644 --- a/applications/solvers/multiphase/multiphaseInterFoam/alphaCourantNo.H +++ b/applications/solvers/multiphase/multiphaseInterFoam/alphaCourantNo.H @@ -39,9 +39,11 @@ scalar meanAlphaCoNum = 0.0; if (mesh.nInternalFaces()) { - scalarField sumPhi = + scalarField sumPhi + ( mixture.nearInterface()().internalField() - *fvc::surfaceSum(mag(phi))().internalField(); + * fvc::surfaceSum(mag(phi))().internalField() + ); alphaCoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue(); diff --git a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.C b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.C index 6256b1cebf..9be378e5c5 100644 --- a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.C +++ b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/alphaContactAngle/alphaContactAngleFvPatchScalarField.C @@ -133,7 +133,11 @@ void alphaContactAngleFvPatchScalarField::write(Ostream& os) const // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -makePatchTypeField(fvPatchScalarField, alphaContactAngleFvPatchScalarField); +makeNonTemplatedPatchTypeField +( + fvPatchScalarField, + alphaContactAngleFvPatchScalarField +); // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C index 22c45bad71..d2bb7d05b4 100644 --- a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C +++ b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C @@ -273,7 +273,7 @@ void Foam::multiphaseMixture::solve() if (nAlphaSubCycles > 1) { - surfaceScalarField rhoPhiSum = 0.0*rhoPhi_; + surfaceScalarField rhoPhiSum(0.0*rhoPhi_); dimensionedScalar totalDeltaT = runTime.deltaT(); for @@ -314,9 +314,11 @@ Foam::tmp Foam::multiphaseMixture::nHatfv surfaceVectorField gradAlphaf = fvc::interpolate(gradAlpha); */ - surfaceVectorField gradAlphaf = + surfaceVectorField gradAlphaf + ( fvc::interpolate(alpha2)*fvc::interpolate(fvc::grad(alpha1)) - - fvc::interpolate(alpha1)*fvc::interpolate(fvc::grad(alpha2)); + - fvc::interpolate(alpha1)*fvc::interpolate(fvc::grad(alpha2)) + ); // Face unit interface normal return gradAlphaf/(mag(gradAlphaf) + deltaN_); @@ -361,9 +363,11 @@ void Foam::multiphaseMixture::correctContactAngle vectorField& nHatPatch = nHatb[patchi]; - vectorField AfHatPatch = + vectorField AfHatPatch + ( mesh_.Sf().boundaryField()[patchi] - /mesh_.magSf().boundaryField()[patchi]; + /mesh_.magSf().boundaryField()[patchi] + ); alphaContactAngleFvPatchScalarField::thetaPropsTable:: const_iterator tp = @@ -396,21 +400,25 @@ void Foam::multiphaseMixture::correctContactAngle scalar thetaR = convertToRad*tp().thetaR(matched); // Calculated the component of the velocity parallel to the wall - vectorField Uwall = + vectorField Uwall + ( U_.boundaryField()[patchi].patchInternalField() - - U_.boundaryField()[patchi]; + - U_.boundaryField()[patchi] + ); Uwall -= (AfHatPatch & Uwall)*AfHatPatch; // Find the direction of the interface parallel to the wall - vectorField nWall = - nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch; + vectorField nWall + ( + nHatPatch - (AfHatPatch & nHatPatch)*AfHatPatch + ); // Normalise nWall nWall /= (mag(nWall) + SMALL); // Calculate Uwall resolved normal to the interface parallel to // the interface - scalarField uwall = nWall & Uwall; + scalarField uwall(nWall & Uwall); theta += (thetaA - thetaR)*tanh(uwall/uTheta); } @@ -418,9 +426,9 @@ void Foam::multiphaseMixture::correctContactAngle // Reset nHatPatch to correspond to the contact angle - scalarField a12 = nHatPatch & AfHatPatch; + scalarField a12(nHatPatch & AfHatPatch); - scalarField b1 = cos(theta); + scalarField b1(cos(theta)); scalarField b2(nHatPatch.size()); @@ -429,10 +437,10 @@ void Foam::multiphaseMixture::correctContactAngle b2[facei] = cos(acos(a12[facei]) - theta[facei]); } - scalarField det = 1.0 - a12*a12; + scalarField det(1.0 - a12*a12); - scalarField a = (b1 - a12*b2)/det; - scalarField b = (b2 - a12*b1)/det; + scalarField a((b1 - a12*b2)/det); + scalarField b((b2 - a12*b1)/det); nHatPatch = a*AfHatPatch + b*nHatPatch; @@ -508,7 +516,7 @@ void Foam::multiphaseMixture::solveAlphas ) ); - surfaceScalarField phic = mag(phi_/mesh_.magSf()); + surfaceScalarField phic(mag(phi_/mesh_.magSf())); phic = min(cAlpha*phic, max(phic)); for (int gCorr=0; gCorr tgradU = fvc::grad(U); - volScalarField G = 2*mut*(tgradU() && dev(symm(tgradU()))); + volScalarField G(2*mut*(tgradU() && dev(symm(tgradU())))); tgradU.clear(); - volScalarField Gcoef = - Cmu*k/sigmak*(g & fvc::grad(rho))/(epsilon + epsilonMin); + volScalarField Gcoef + ( + Cmu*k/sigmak*(g & fvc::grad(rho))/(epsilon + epsilonMin) + ); #include "wallFunctions.H" diff --git a/applications/solvers/multiphase/settlingFoam/pEqn.H b/applications/solvers/multiphase/settlingFoam/pEqn.H index 806909ffe7..f359e95254 100644 --- a/applications/solvers/multiphase/settlingFoam/pEqn.H +++ b/applications/solvers/multiphase/settlingFoam/pEqn.H @@ -1,4 +1,4 @@ -volScalarField rAU = 1.0/UEqn.A(); +volScalarField rAU(1.0/UEqn.A()); surfaceScalarField rAUf ( diff --git a/applications/solvers/multiphase/settlingFoam/plasticViscosity.H b/applications/solvers/multiphase/settlingFoam/plasticViscosity.H index cbbeb7f3f5..2797205a75 100644 --- a/applications/solvers/multiphase/settlingFoam/plasticViscosity.H +++ b/applications/solvers/multiphase/settlingFoam/plasticViscosity.H @@ -5,11 +5,13 @@ volScalarField plasticViscosity const volScalarField& alpha ) { - return + tmp tfld ( plasticViscosityCoeff* ( pow(10.0, plasticViscosityExponent*alpha + SMALL) - scalar(1) ) ); + + return tfld(); } diff --git a/applications/solvers/multiphase/settlingFoam/wallFunctions.H b/applications/solvers/multiphase/settlingFoam/wallFunctions.H index decccfbff5..952d3c4591 100644 --- a/applications/solvers/multiphase/settlingFoam/wallFunctions.H +++ b/applications/solvers/multiphase/settlingFoam/wallFunctions.H @@ -34,11 +34,13 @@ { const scalarField& rhow = rho.boundaryField()[patchi]; - const scalarField muw = mul.boundaryField()[patchi]; + const scalarField muw(mul.boundaryField()[patchi]); const scalarField& mutw = mut.boundaryField()[patchi]; - scalarField magFaceGradU = - mag(U.boundaryField()[patchi].snGrad()); + scalarField magFaceGradU + ( + mag(U.boundaryField()[patchi].snGrad()) + ); forAll(curPatch, facei) { diff --git a/applications/solvers/multiphase/settlingFoam/wallViscosity.H b/applications/solvers/multiphase/settlingFoam/wallViscosity.H index 8f64c27b3c..9516bb0fea 100644 --- a/applications/solvers/multiphase/settlingFoam/wallViscosity.H +++ b/applications/solvers/multiphase/settlingFoam/wallViscosity.H @@ -13,7 +13,7 @@ { const scalarField& rhow = rho.boundaryField()[patchi]; - const scalarField muw = mul.boundaryField()[patchi]; + const scalarField muw(mul.boundaryField()[patchi]); scalarField& mutw = mut.boundaryField()[patchi]; forAll(curPatch, facei) diff --git a/applications/solvers/multiphase/settlingFoam/yieldStress.H b/applications/solvers/multiphase/settlingFoam/yieldStress.H index 6b1b734f50..cb0415c66f 100644 --- a/applications/solvers/multiphase/settlingFoam/yieldStress.H +++ b/applications/solvers/multiphase/settlingFoam/yieldStress.H @@ -6,7 +6,7 @@ volScalarField yieldStress const volScalarField& alpha ) { - return + tmp tfld ( yieldStressCoeff* ( @@ -14,4 +14,6 @@ volScalarField yieldStress - pow(10.0, yieldStressExponent*yieldStressOffset) ) ); + + return tfld(); } diff --git a/applications/solvers/multiphase/twoLiquidMixingFoam/pEqn.H b/applications/solvers/multiphase/twoLiquidMixingFoam/pEqn.H index ca11063449..08f00fc185 100644 --- a/applications/solvers/multiphase/twoLiquidMixingFoam/pEqn.H +++ b/applications/solvers/multiphase/twoLiquidMixingFoam/pEqn.H @@ -1,6 +1,6 @@ { - volScalarField rAU = 1.0/UEqn.A(); - surfaceScalarField rAUf = fvc::interpolate(rAU); + volScalarField rAU(1.0/UEqn.A()); + surfaceScalarField rAUf(fvc::interpolate(rAU)); U = rAU*UEqn.H(); surfaceScalarField phiU diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H b/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H index 6eb2fe8850..b07f3e70dd 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H @@ -3,7 +3,7 @@ fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime); { { - volTensorField gradUaT = fvc::grad(Ua)().T(); + volTensorField gradUaT(fvc::grad(Ua)().T()); if (kineticTheory.on()) { @@ -26,9 +26,11 @@ fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime); Rca -= ((kineticTheory.lambda()/rhoa)*tr(gradUaT))*tensor(I); } - surfaceScalarField phiRa = + surfaceScalarField phiRa + ( -fvc::interpolate(nuEffa)*mesh.magSf()*fvc::snGrad(alpha) - /fvc::interpolate(alpha + scalar(0.001)); + /fvc::interpolate(alpha + scalar(0.001)) + ); UaEqn = ( @@ -56,16 +58,18 @@ fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime); } { - volTensorField gradUbT = fvc::grad(Ub)().T(); + volTensorField gradUbT(fvc::grad(Ub)().T()); volTensorField Rcb ( "Rcb", ((2.0/3.0)*I)*(k + nuEffb*tr(gradUbT)) - nuEffb*gradUbT ); - surfaceScalarField phiRb = + surfaceScalarField phiRb + ( -fvc::interpolate(nuEffb)*mesh.magSf()*fvc::snGrad(beta) - /fvc::interpolate(beta + scalar(0.001)); + /fvc::interpolate(beta + scalar(0.001)) + ); UbEqn = ( diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/alphaEqn.H b/applications/solvers/multiphase/twoPhaseEulerFoam/alphaEqn.H index 47057c0efb..cac8e52a8e 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/alphaEqn.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/alphaEqn.H @@ -2,13 +2,13 @@ word scheme("div(phi,alpha)"); word schemer("div(phir,alpha)"); - surfaceScalarField phic = phi; - surfaceScalarField phir = phia - phib; + surfaceScalarField phic(phi); + surfaceScalarField phir(phia - phib); if (g0.value() > 0.0) { - surfaceScalarField alphaf = fvc::interpolate(alpha); - surfaceScalarField phipp = ppMagf*fvc::snGrad(alpha)*mesh.magSf(); + surfaceScalarField alphaf(fvc::interpolate(alpha)); + surfaceScalarField phipp(ppMagf*fvc::snGrad(alpha)*mesh.magSf()); phir += phipp; phic += fvc::interpolate(alpha)*phipp; } diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/createFields.H b/applications/solvers/multiphase/twoPhaseEulerFoam/createFields.H index e422fe845b..00581cd302 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/createFields.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/createFields.H @@ -132,15 +132,19 @@ Info<< "Calculating field DDtUa and DDtUb\n" << endl; - volVectorField DDtUa = + volVectorField DDtUa + ( fvc::ddt(Ua) + fvc::div(phia, Ua) - - fvc::div(phia)*Ua; + - fvc::div(phia)*Ua + ); - volVectorField DDtUb = + volVectorField DDtUb + ( fvc::ddt(Ub) + fvc::div(phib, Ub) - - fvc::div(phib)*Ub; + - fvc::div(phib)*Ub + ); Info<< "Calculating field g.h\n" << endl; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C index da4c723a32..2b68640f23 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C @@ -68,7 +68,7 @@ Foam::tmp Foam::Ergun::K const volScalarField& Ur ) const { - volScalarField beta = max(scalar(1) - alpha_, scalar(1.0e-6)); + volScalarField beta(max(scalar(1) - alpha_, scalar(1.0e-6))); return 150.0*alpha_*phaseb_.nu()*phaseb_.rho() diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C index e728d5d967..4ff7d0fc13 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C @@ -68,9 +68,9 @@ Foam::tmp Foam::Gibilaro::K const volScalarField& Ur ) const { - volScalarField beta = max(scalar(1) - alpha_, scalar(1.0e-6)); - volScalarField bp = pow(beta, -2.8); - volScalarField Re = max(beta*Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3)); + volScalarField beta(max(scalar(1) - alpha_, scalar(1.0e-6))); + volScalarField bp(pow(beta, -2.8)); + volScalarField Re(max(beta*Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3))); return (17.3/Re + scalar(0.336))*phaseb_.rho()*Ur*bp/phasea_.d(); } diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C index da72226edb..a99f975330 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C @@ -68,12 +68,12 @@ Foam::tmp Foam::GidaspowErgunWenYu::K const volScalarField& Ur ) const { - volScalarField beta = max(scalar(1) - alpha_, scalar(1.0e-6)); + volScalarField beta(max(scalar(1) - alpha_, scalar(1.0e-6))); - volScalarField bp = pow(beta, -2.65); - volScalarField Re = max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3)); + volScalarField bp(pow(beta, -2.65)); + volScalarField Re(max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3))); - volScalarField Cds = 24.0*(1.0 + 0.15*pow(Re, 0.687))/Re; + volScalarField Cds(24.0*(1.0 + 0.15*pow(Re, 0.687))/Re); forAll(Re, celli) { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C index 5f9f2054db..1da8e0c547 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C @@ -68,11 +68,11 @@ Foam::tmp Foam::GidaspowSchillerNaumann::K const volScalarField& Ur ) const { - volScalarField beta = max(scalar(1) - alpha_, scalar(1e-6)); - volScalarField bp = pow(beta, -2.65); + volScalarField beta(max(scalar(1) - alpha_, scalar(1e-6))); + volScalarField bp(pow(beta, -2.65)); - volScalarField Re = max(beta*Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3)); - volScalarField Cds = 24.0*(scalar(1) + 0.15*pow(Re, 0.687))/Re; + volScalarField Re(max(beta*Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3))); + volScalarField Cds(24.0*(scalar(1) + 0.15*pow(Re, 0.687))/Re); forAll(Re, celli) { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C index 1e99ae05f5..cdfaac389e 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C @@ -68,8 +68,8 @@ Foam::tmp Foam::SchillerNaumann::K const volScalarField& Ur ) const { - volScalarField Re = max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3)); - volScalarField Cds = 24.0*(scalar(1) + 0.15*pow(Re, 0.687))/Re; + volScalarField Re(max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3))); + volScalarField Cds(24.0*(scalar(1) + 0.15*pow(Re, 0.687))/Re); forAll(Re, celli) { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C index 517b21cfe7..45989dbf3f 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C @@ -68,9 +68,9 @@ Foam::tmp Foam::SyamlalOBrien::K const volScalarField& Ur ) const { - volScalarField beta = max(scalar(1) - alpha_, scalar(1.0e-6)); - volScalarField A = pow(beta, 4.14); - volScalarField B = 0.8*pow(beta, 1.28); + volScalarField beta(max(scalar(1) - alpha_, scalar(1.0e-6))); + volScalarField A(pow(beta, 4.14)); + volScalarField B(0.8*pow(beta, 1.28)); forAll (beta, celli) { @@ -80,14 +80,17 @@ Foam::tmp Foam::SyamlalOBrien::K } } - volScalarField Re = max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3)); + volScalarField Re(max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3))); - volScalarField Vr = 0.5* + volScalarField Vr ( - A - 0.06*Re + sqrt(sqr(0.06*Re) + 0.12*Re*(2.0*B - A) + sqr(A)) + 0.5* + ( + A - 0.06*Re + sqrt(sqr(0.06*Re) + 0.12*Re*(2.0*B - A) + sqr(A)) + ) ); - volScalarField Cds = sqr(0.63 + 4.8*sqrt(Vr/Re)); + volScalarField Cds(sqr(0.63 + 4.8*sqrt(Vr/Re))); return 0.75*Cds*phaseb_.rho()*Ur/(phasea_.d()*sqr(Vr)); } diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C index 164d82a4a5..bfd5eaadbd 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C @@ -68,11 +68,11 @@ Foam::tmp Foam::WenYu::K const volScalarField& Ur ) const { - volScalarField beta = max(scalar(1) - alpha_, scalar(1.0e-6)); - volScalarField bp = pow(beta, -2.65); + volScalarField beta(max(scalar(1) - alpha_, scalar(1.0e-6))); + volScalarField bp(pow(beta, -2.65)); - volScalarField Re = max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3)); - volScalarField Cds = 24.0*(scalar(1) + 0.15*pow(Re, 0.687))/Re; + volScalarField Re(max(Ur*phasea_.d()/phaseb_.nu(), scalar(1.0e-3))); + volScalarField Cds(24.0*(scalar(1) + 0.15*pow(Re, 0.687))/Re); forAll(Re, celli) { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kEpsilon.H b/applications/solvers/multiphase/twoPhaseEulerFoam/kEpsilon.H index 3d53ef8684..8f85e12052 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kEpsilon.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kEpsilon.H @@ -6,7 +6,7 @@ if (turbulence) } tmp tgradUb = fvc::grad(Ub); - volScalarField G = 2*nutb*(tgradUb() && dev(symm(tgradUb()))); + volScalarField G(2*nutb*(tgradUb() && dev(symm(tgradUb())))); tgradUb.clear(); #include "wallFunctions.H" diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/conductivityModel/HrenyaSinclair/HrenyaSinclairConductivity.C b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/conductivityModel/HrenyaSinclair/HrenyaSinclairConductivity.C index c1eb55e14e..d95a2d2831 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/conductivityModel/HrenyaSinclair/HrenyaSinclairConductivity.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/conductivityModel/HrenyaSinclair/HrenyaSinclairConductivity.C @@ -75,8 +75,10 @@ Foam::tmp Foam::HrenyaSinclairConductivity::kappa { const scalar sqrtPi = sqrt(constant::mathematical::pi); - volScalarField lamda = - scalar(1) + da/(6.0*sqrt(2.0)*(alpha + scalar(1.0e-5)))/L_; + volScalarField lamda + ( + scalar(1) + da/(6.0*sqrt(2.0)*(alpha + scalar(1.0e-5)))/L_ + ); return rhoa*da*sqrt(Theta)* ( diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C index ee2300d91e..cb03cb2761 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C @@ -202,16 +202,16 @@ void Foam::kineticTheoryModel::solve(const volTensorField& gradUat) const scalar sqrtPi = sqrt(constant::mathematical::pi); - surfaceScalarField phi = 1.5*rhoa_*phia_*fvc::interpolate(alpha_); + surfaceScalarField phi(1.5*rhoa_*phia_*fvc::interpolate(alpha_)); - volTensorField dU = gradUat.T();//fvc::grad(Ua_); - volSymmTensorField D = symm(dU); + volTensorField dU(gradUat.T()); //fvc::grad(Ua_); + volSymmTensorField D(symm(dU)); // NB, drag = K*alpha*beta, // (the alpha and beta has been extracted from the drag function for // numerical reasons) - volScalarField Ur = mag(Ua_ - Ub_); - volScalarField betaPrim = alpha_*(1.0 - alpha_)*draga_.K(Ur); + volScalarField Ur(mag(Ua_ - Ub_)); + volScalarField betaPrim(alpha_*(1.0 - alpha_)*draga_.K(Ur)); // Calculating the radial distribution function (solid volume fraction is // limited close to the packing limit, but this needs improvements) @@ -223,12 +223,15 @@ void Foam::kineticTheoryModel::solve(const volTensorField& gradUat) ); // particle pressure - coefficient in front of Theta (Eq. 3.22, p. 45) - volScalarField PsCoeff = granularPressureModel_->granularPressureCoeff + volScalarField PsCoeff ( - alpha_, - gs0_, - rhoa_, - e_ + granularPressureModel_->granularPressureCoeff + ( + alpha_, + gs0_, + rhoa_, + e_ + ) ); // 'thermal' conductivity (Table 3.3, p. 49) @@ -245,23 +248,27 @@ void Foam::kineticTheoryModel::solve(const volTensorField& gradUat) ); dimensionedScalar TsmallSqrt = sqrt(Tsmall); - volScalarField ThetaSqrt = sqrt(Theta_); + volScalarField ThetaSqrt(sqrt(Theta_)); // dissipation (Eq. 3.24, p.50) - volScalarField gammaCoeff = - 12.0*(1.0 - sqr(e_))*sqr(alpha_)*rhoa_*gs0_*(1.0/da_)*ThetaSqrt/sqrtPi; + volScalarField gammaCoeff + ( + 12.0*(1.0 - sqr(e_))*sqr(alpha_)*rhoa_*gs0_*(1.0/da_)*ThetaSqrt/sqrtPi + ); // Eq. 3.25, p. 50 Js = J1 - J2 - volScalarField J1 = 3.0*betaPrim; - volScalarField J2 = + volScalarField J1(3.0*betaPrim); + volScalarField J2 + ( 0.25*sqr(betaPrim)*da_*sqr(Ur) - /(max(alpha_, 1e-6)*rhoa_*sqrtPi*(ThetaSqrt + TsmallSqrt)); + / (max(alpha_, 1e-6)*rhoa_*sqrtPi*(ThetaSqrt + TsmallSqrt)) + ); // bulk viscosity p. 45 (Lun et al. 1984). lambda_ = (4.0/3.0)*sqr(alpha_)*rhoa_*da_*gs0_*(1.0+e_)*ThetaSqrt/sqrtPi; // stress tensor, Definitions, Table 3.1, p. 43 - volSymmTensorField tau = 2.0*mua_*D + (lambda_ - (2.0/3.0)*mua_)*tr(D)*I; + volSymmTensorField tau(2.0*mua_*D + (lambda_ - (2.0/3.0)*mua_)*tr(D)*I); if (!equilibrium_) { @@ -289,27 +296,32 @@ void Foam::kineticTheoryModel::solve(const volTensorField& gradUat) { // equilibrium => dissipation == production // Eq. 4.14, p.82 - volScalarField K1 = 2.0*(1.0 + e_)*rhoa_*gs0_; - volScalarField K3 = 0.5*da_*rhoa_* + volScalarField K1(2.0*(1.0 + e_)*rhoa_*gs0_); + volScalarField K3 + ( + 0.5*da_*rhoa_* ( (sqrtPi/(3.0*(3.0-e_))) *(1.0 + 0.4*(1.0 + e_)*(3.0*e_ - 1.0)*alpha_*gs0_) +1.6*alpha_*gs0_*(1.0 + e_)/sqrtPi - ); + ) + ); - volScalarField K2 = - 4.0*da_*rhoa_*(1.0 + e_)*alpha_*gs0_/(3.0*sqrtPi) - 2.0*K3/3.0; + volScalarField K2 + ( + 4.0*da_*rhoa_*(1.0 + e_)*alpha_*gs0_/(3.0*sqrtPi) - 2.0*K3/3.0 + ); - volScalarField K4 = 12.0*(1.0 - sqr(e_))*rhoa_*gs0_/(da_*sqrtPi); + volScalarField K4(12.0*(1.0 - sqr(e_))*rhoa_*gs0_/(da_*sqrtPi)); - volScalarField trD = tr(D); - volScalarField tr2D = sqr(trD); - volScalarField trD2 = tr(D & D); + volScalarField trD(tr(D)); + volScalarField tr2D(sqr(trD)); + volScalarField trD2(tr(D & D)); - volScalarField t1 = K1*alpha_ + rhoa_; - volScalarField l1 = -t1*trD; - volScalarField l2 = sqr(t1)*tr2D; - volScalarField l3 = 4.0*K4*max(alpha_, 1e-6)*(2.0*K3*trD2 + K2*tr2D); + volScalarField t1(K1*alpha_ + rhoa_); + volScalarField l1(-t1*trD); + volScalarField l2(sqr(t1)*tr2D); + volScalarField l3(4.0*K4*max(alpha_, 1e-6)*(2.0*K3*trD2 + K2*tr2D)); Theta_ = sqr((l1 + sqrt(l2 + l3))/(2.0*(alpha_ + 1.0e-4)*K4)); } @@ -317,14 +329,17 @@ void Foam::kineticTheoryModel::solve(const volTensorField& gradUat) Theta_.max(1.0e-15); Theta_.min(1.0e+3); - volScalarField pf = frictionalStressModel_->frictionalPressure + volScalarField pf ( - alpha_, - alphaMinFriction_, - alphaMax_, - Fr_, - eta_, - p_ + frictionalStressModel_->frictionalPressure + ( + alpha_, + alphaMinFriction_, + alphaMax_, + Fr_, + eta_, + p_ + ) ); PsCoeff += pf/(Theta_+Tsmall); @@ -336,23 +351,26 @@ void Foam::kineticTheoryModel::solve(const volTensorField& gradUat) pa_ = PsCoeff*Theta_; // frictional shear stress, Eq. 3.30, p. 52 - volScalarField muf = frictionalStressModel_->muf + volScalarField muf ( - alpha_, - alphaMax_, - pf, - D, - phi_ + frictionalStressModel_->muf + ( + alpha_, + alphaMax_, + pf, + D, + phi_ + ) ); - // add frictional stress + // add frictional stress mua_ += muf; mua_.min(1.0e+2); mua_.max(0.0); Info<< "kinTheory: max(Theta) = " << max(Theta_).value() << endl; - volScalarField ktn = mua_/rhoa_; + volScalarField ktn(mua_/rhoa_); Info<< "kinTheory: min(nua) = " << min(ktn).value() << ", max(nua) = " << max(ktn).value() << endl; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/viscosityModel/HrenyaSinclair/HrenyaSinclairViscosity.C b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/viscosityModel/HrenyaSinclair/HrenyaSinclairViscosity.C index e52a3b204e..620ca79b69 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/viscosityModel/HrenyaSinclair/HrenyaSinclairViscosity.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/viscosityModel/HrenyaSinclair/HrenyaSinclairViscosity.C @@ -78,8 +78,10 @@ Foam::kineticTheoryModels::HrenyaSinclairViscosity::mua { const scalar sqrtPi = sqrt(constant::mathematical::pi); - volScalarField lamda = - scalar(1) + da/(6.0*sqrt(2.0)*(alpha + scalar(1.0e-5)))/L_; + volScalarField lamda + ( + scalar(1) + da/(6.0*sqrt(2.0)*(alpha + scalar(1.0e-5)))/L_ + ); return rhoa*da*sqrt(Theta)* ( diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/liftDragCoeffs.H b/applications/solvers/multiphase/twoPhaseEulerFoam/liftDragCoeffs.H index 5e0caf1f09..b614a4dfbe 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/liftDragCoeffs.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/liftDragCoeffs.H @@ -1,18 +1,18 @@ - volVectorField Ur = Ua - Ub; - volScalarField magUr = mag(Ur); + volVectorField Ur(Ua - Ub); + volScalarField magUr(mag(Ur)); - volScalarField Ka = draga->K(magUr); - volScalarField K = Ka; + volScalarField Ka(draga->K(magUr)); + volScalarField K(Ka); if (dragPhase == "b") { - volScalarField Kb = dragb->K(magUr); + volScalarField Kb(dragb->K(magUr)); K = Kb; } else if (dragPhase == "blended") { - volScalarField Kb = dragb->K(magUr); + volScalarField Kb(dragb->K(magUr)); K = (beta*Ka + alpha*Kb); } - volVectorField liftCoeff = Cl*(beta*rhob + alpha*rhoa)*(Ur ^ fvc::curl(U)); + volVectorField liftCoeff(Cl*(beta*rhob + alpha*rhoa)*(Ur ^ fvc::curl(U))); diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H b/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H index f92944a414..82fc1280c1 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H @@ -1,21 +1,23 @@ { - surfaceScalarField alphaf = fvc::interpolate(alpha); - surfaceScalarField betaf = scalar(1) - alphaf; + surfaceScalarField alphaf(fvc::interpolate(alpha)); + surfaceScalarField betaf(scalar(1) - alphaf); - volScalarField rUaA = 1.0/UaEqn.A(); - volScalarField rUbA = 1.0/UbEqn.A(); + volScalarField rUaA(1.0/UaEqn.A()); + volScalarField rUbA(1.0/UbEqn.A()); phia == (fvc::interpolate(Ua) & mesh.Sf()); phib == (fvc::interpolate(Ub) & mesh.Sf()); rUaAf = fvc::interpolate(rUaA); - surfaceScalarField rUbAf = fvc::interpolate(rUbA); + surfaceScalarField rUbAf(fvc::interpolate(rUbA)); Ua = rUaA*UaEqn.H(); Ub = rUbA*UbEqn.H(); - surfaceScalarField phiDraga = - fvc::interpolate(beta/rhoa*K*rUaA)*phib + rUaAf*(g & mesh.Sf()); + surfaceScalarField phiDraga + ( + fvc::interpolate(beta/rhoa*K*rUaA)*phib + rUaAf*(g & mesh.Sf()) + ); if (g0.value() > 0.0) { @@ -27,8 +29,10 @@ phiDraga -= rUaAf*fvc::snGrad(kineticTheory.pa()/rhoa)*mesh.magSf(); } - surfaceScalarField phiDragb = - fvc::interpolate(alpha/rhob*K*rUbA)*phia + rUbAf*(g & mesh.Sf()); + surfaceScalarField phiDragb + ( + fvc::interpolate(alpha/rhob*K*rUbA)*phia + rUbAf*(g & mesh.Sf()) + ); // Fix for gravity on outlet boundary. forAll(p.boundaryField(), patchi) @@ -66,7 +70,7 @@ if (nonOrth == nNonOrthCorr) { - surfaceScalarField SfGradp = pEqn.flux()/Dp; + surfaceScalarField SfGradp(pEqn.flux()/Dp); phia -= rUaAf*SfGradp/rhoa; phib -= rUbAf*SfGradp/rhob; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/packingLimiter.H b/applications/solvers/multiphase/twoPhaseEulerFoam/packingLimiter.H index ba6dcf7ce6..e5d9d86a55 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/packingLimiter.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/packingLimiter.H @@ -1,11 +1,11 @@ if (packingLimiter) { // Calculating exceeding volume fractions - volScalarField alphaEx = max(alpha - alphaMax, scalar(0)); + volScalarField alphaEx(max(alpha - alphaMax, scalar(0))); // Finding neighbouring cells of the whole domain labelListList neighbour = mesh.cellCells(); - scalarField cellVolumes = mesh.cellVolumes(); + scalarField cellVolumes(mesh.cellVolumes()); forAll (alphaEx, celli) {