diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H b/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H index 9f0417045b..c32855904a 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/UEqns.H @@ -4,12 +4,24 @@ fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime); { { volTensorField gradUaT = fvc::grad(Ua)().T(); + + if (kineticTheory.on()) + { + kineticTheory.solve(gradUaT); + nuEffa = kineticTheory.mua()/rhoa; + } + else // If not using kinetic theory is using Ct model + { + nuEffa = sqr(Ct)*nutb + nua; + } + volTensorField Rca ( "Rca", ((2.0/3.0)*I)*(sqr(Ct)*k + nuEffa*tr(gradUaT)) - nuEffa*gradUaT ); + if (kineticTheory.on()) { Rca -= ((kineticTheory.lambda()/rhoa)*tr(gradUaT))*tensor(I); diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kEpsilon.H b/applications/solvers/multiphase/twoPhaseEulerFoam/kEpsilon.H new file mode 100644 index 0000000000..dab75ee07d --- /dev/null +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kEpsilon.H @@ -0,0 +1,62 @@ +if(turbulence) +{ + if (mesh.changing()) + { + y.correct(); + } + + tmp tgradUb = fvc::grad(Ub); + volScalarField G = 2*nutb*(tgradUb() && dev(symm(tgradUb()))); + tgradUb.clear(); + + #include "wallFunctions.H" + + // Dissipation equation + fvScalarMatrix epsEqn + ( + fvm::ddt(beta, epsilon) + + fvm::div(phib, epsilon) + - fvm::laplacian + ( + alphaEps*nuEffb, epsilon, + "laplacian(DepsilonEff,epsilon)" + ) + == + C1*beta*G*epsilon/k + - fvm::Sp(C2*beta*epsilon/k, epsilon) + ); + + #include "wallDissipation.H" + + epsEqn.relax(); + epsEqn.solve(); + + epsilon.max(dimensionedScalar("zero", epsilon.dimensions(), 1.0e-15)); + + + // Turbulent kinetic energy equation + fvScalarMatrix kEqn + ( + fvm::ddt(beta, k) + + fvm::div(phib, k) + - fvm::laplacian + ( + alphak*nuEffb, k, + "laplacian(DkEff,k)" + ) + == + beta*G + - fvm::Sp(beta*epsilon/k, k) + ); + kEqn.relax(); + kEqn.solve(); + + k.max(dimensionedScalar("zero", k.dimensions(), 1.0e-8)); + + //- Re-calculate turbulence viscosity + nutb = Cmu*sqr(k)/epsilon; + + #include "wallViscosity.H" +} + +nuEffb = nutb + nub; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.C b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.C index af326ffd5b..bcd41ac268 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.C @@ -102,7 +102,7 @@ Foam::tmp Foam::JohnsonJacksonFrictionalStress::muf const volScalarField& alpha, const dimensionedScalar& alphaMax, const volScalarField& pf, - const volTensorField& D, + const volSymmTensorField& D, const dimensionedScalar& phi ) const { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.H b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.H index 08733247c5..d758f2f85c 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/JohnsonJackson/JohnsonJacksonFrictionalStress.H @@ -93,7 +93,7 @@ public: const volScalarField& alpha, const dimensionedScalar& alphaMax, const volScalarField& pf, - const volTensorField& D, + const volSymmTensorField& D, const dimensionedScalar& phi ) const; }; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C index dbc658f1d4..0f39146345 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.C @@ -99,7 +99,7 @@ Foam::tmp Foam::SchaefferFrictionalStress::muf const volScalarField& alpha, const dimensionedScalar& alphaMax, const volScalarField& pf, - const volTensorField& D, + const volSymmTensorField& D, const dimensionedScalar& phi ) const { @@ -124,9 +124,9 @@ Foam::tmp Foam::SchaefferFrictionalStress::muf volScalarField& muff = tmuf(); - forAll(D, celli) + forAll (D, celli) { - if (alpha[celli] > alphaMax.value()-5e-2) + if (alpha[celli] > alphaMax.value() - 5e-2) { muff[celli] = 0.5*pf[celli]*sin(phi.value()) diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.H b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.H index e9b0c4efd7..89b0e25ded 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/Schaeffer/SchaefferFrictionalStress.H @@ -93,7 +93,7 @@ public: const volScalarField& alpha, const dimensionedScalar& alphaMax, const volScalarField& pf, - const volTensorField& D, + const volSymmTensorField& D, const dimensionedScalar& phi ) const; }; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/frictionalStressModel/frictionalStressModel.H b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/frictionalStressModel/frictionalStressModel.H index 8176e6aac2..b496ee9cf6 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/frictionalStressModel/frictionalStressModel.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/frictionalStressModel/frictionalStressModel/frictionalStressModel.H @@ -48,7 +48,7 @@ namespace Foam class frictionalStressModel { - // Private Member Functions + // Private member functions //- Disallow default bitwise copy construct frictionalStressModel(const frictionalStressModel&); @@ -127,7 +127,7 @@ public: const volScalarField& alpha, const dimensionedScalar& alphaMax, const volScalarField& pf, - const volTensorField& D, + const volSymmTensorField& D, const dimensionedScalar& phi ) const = 0; }; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C index abdd6080be..ee2300d91e 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C @@ -56,12 +56,13 @@ Foam::kineticTheoryModel::kineticTheoryModel "kineticTheoryProperties", Ua_.time().constant(), Ua_.mesh(), - IOobject::MUST_READ_IF_MODIFIED, + IOobject::MUST_READ, IOobject::NO_WRITE ) ), kineticTheory_(kineticTheoryProperties_.lookup("kineticTheory")), equilibrium_(kineticTheoryProperties_.lookup("equilibrium")), + viscosityModel_ ( kineticTheoryModels::viscosityModel::New @@ -192,24 +193,19 @@ Foam::kineticTheoryModel::~kineticTheoryModel() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -void Foam::kineticTheoryModel::solve() +void Foam::kineticTheoryModel::solve(const volTensorField& gradUat) { if (!kineticTheory_) { return; } - word scheme("div(phi,Theta)"); - - volScalarField alpha = alpha_; - alpha.max(1.0e-6); const scalar sqrtPi = sqrt(constant::mathematical::pi); surfaceScalarField phi = 1.5*rhoa_*phia_*fvc::interpolate(alpha_); - volTensorField dU = fvc::grad(Ua_); - volTensorField dUT = dU.T(); - volTensorField D = 0.5*(dU + dUT); + 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 @@ -220,45 +216,52 @@ void Foam::kineticTheoryModel::solve() // Calculating the radial distribution function (solid volume fraction is // limited close to the packing limit, but this needs improvements) // The solution is higly unstable close to the packing limit. - gs0_ = radialModel_->g0(min(alpha, alphaMax_-1.0e-2), alphaMax_); + gs0_ = radialModel_->g0 + ( + min(max(alpha_, 1e-6), alphaMax_ - 0.01), + alphaMax_ + ); // particle pressure - coefficient in front of Theta (Eq. 3.22, p. 45) - volScalarField PsCoeff = - granularPressureModel_->granularPressureCoeff(alpha_,gs0_,rhoa_,e_ ); + volScalarField PsCoeff = granularPressureModel_->granularPressureCoeff + ( + alpha_, + gs0_, + rhoa_, + e_ + ); + + // 'thermal' conductivity (Table 3.3, p. 49) + kappa_ = conductivityModel_->kappa(alpha_, Theta_, gs0_, rhoa_, da_, e_); + + // particle viscosity (Table 3.2, p.47) + mua_ = viscosityModel_->mua(alpha_, Theta_, gs0_, rhoa_, da_, e_); dimensionedScalar Tsmall ( "small", - dimensionSet(0,2,-2,0,0,0,0), + dimensionSet(0 , 2 ,-2 ,0 , 0, 0, 0), 1.0e-6 ); dimensionedScalar TsmallSqrt = sqrt(Tsmall); volScalarField ThetaSqrt = sqrt(Theta_); - // 'thermal' conductivity (Table 3.3, p. 49) - kappa_ = conductivityModel_->kappa(alpha, Theta_, gs0_, rhoa_, da_, e_); - - // particle viscosity (Table 3.2, p.47) - mua_ = viscosityModel_->mua(alpha, Theta_, gs0_, rhoa_, da_, e_); - // dissipation (Eq. 3.24, p.50) volScalarField gammaCoeff = - 12.0*(1.0 - e_*e_)*sqr(alpha)*rhoa_*gs0_*(1.0/da_) - *ThetaSqrt/sqrtPi; + 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 = 0.25*sqr(betaPrim)*da_*sqr(Ur) - /(alpha*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 - volTensorField 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_) { @@ -268,8 +271,8 @@ void Foam::kineticTheoryModel::solve() // wrong sign infront of laplacian fvScalarMatrix ThetaEqn ( - fvm::ddt(1.5*alpha*rhoa_, Theta_) - + fvm::div(phi, Theta_, scheme) + fvm::ddt(1.5*alpha_*rhoa_, Theta_) + + fvm::div(phi, Theta_, "div(phi,Theta)") == fvm::SuSp(-((PsCoeff*I) && dU), Theta_) + (tau && dU) @@ -290,33 +293,31 @@ void Foam::kineticTheoryModel::solve() 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 + *(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; + 4.0*da_*rhoa_*(1.0 + e_)*alpha_*gs0_/(3.0*sqrtPi) - 2.0*K3/3.0; - volScalarField K4 = 12.0*(1.0 - e_*e_)*rhoa_*gs0_/(da_*sqrtPi); + volScalarField K4 = 12.0*(1.0 - sqr(e_))*rhoa_*gs0_/(da_*sqrtPi); volScalarField trD = tr(D); - volTensorField D2 = D & D; - volScalarField tr2D = trD*trD; - volScalarField trD2 = tr(D2); + volScalarField tr2D = sqr(trD); + volScalarField trD2 = tr(D & D); - volScalarField t1 = K1*alpha + rhoa_; + volScalarField t1 = K1*alpha_ + rhoa_; volScalarField l1 = -t1*trD; volScalarField l2 = sqr(t1)*tr2D; - volScalarField l3 = 4.0*K4*alpha*(2.0*K3*trD2 + K2*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)); + Theta_ = sqr((l1 + sqrt(l2 + l3))/(2.0*(alpha_ + 1.0e-4)*K4)); } Theta_.max(1.0e-15); Theta_.min(1.0e+3); - volScalarField pf = - frictionalStressModel_->frictionalPressure + volScalarField pf = frictionalStressModel_->frictionalPressure ( alpha_, alphaMinFriction_, @@ -344,13 +345,11 @@ void Foam::kineticTheoryModel::solve() phi_ ); - // add frictional stress for alpha > alphaMinFriction - mua_ = viscosityModel_->mua(alpha, Theta_, gs0_, rhoa_, da_, e_) + muf; + // add frictional stress + mua_ += muf; mua_.min(1.0e+2); mua_.max(0.0); - lambda_ = (4.0/3.0)*sqr(alpha_)*rhoa_*da_*gs0_*(1.0 + e_)*ThetaSqrt/sqrtPi; - Info<< "kinTheory: max(Theta) = " << max(Theta_).value() << endl; volScalarField ktn = mua_/rhoa_; diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H index 2659d8c1d6..2863c2d038 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H @@ -156,7 +156,7 @@ public: // Member Functions - void solve(); + void solve(const volTensorField& gradUat); bool on() const { diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H b/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H index b2bc0d817d..39a133d39b 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/pEqn.H @@ -5,6 +5,9 @@ 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); @@ -47,11 +50,10 @@ surfaceScalarField Dp ( - "(rho*(1|A(U)))", - alphaf*rUaAf/rhoa + betaf*rUbAf/rhob + "(rho*(1|A(U)))", alphaf*rUaAf/rhoa + betaf*rUbAf/rhob ); - for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + for(int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) { fvScalarMatrix pEqn ( diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseEulerFoam.C b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseEulerFoam.C index b3e0eb0009..bd0310765d 100644 --- a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseEulerFoam.C +++ b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseEulerFoam.C @@ -92,11 +92,6 @@ int main(int argc, char *argv[]) #include "kEpsilon.H" - if (kineticTheory.on()) - { - kineticTheory.solve(); - nuEffa += kineticTheory.mua()/rhoa; - } #include "write.H" Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" diff --git a/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake b/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake index 1e2ebd6c43..cd3a58ffed 100755 --- a/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake +++ b/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake @@ -7,9 +7,9 @@ then case "$ParaView_VERSION" in 3* | git) - if [ ! -d "${PV_PLUGIN_PATH}" ] + if [ ! -n "${PV_PLUGIN_PATH}" ] then - echo "$0 : PV_PLUGIN_PATH not a valid directory." + echo "$0 : PV_PLUGIN_PATH not a valid." exit 1 fi diff --git a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldFunctions.H b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldFunctions.H index e66d423fa7..1f846416a5 100644 --- a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldFunctions.H +++ b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldFunctions.H @@ -126,8 +126,6 @@ BINARY_FUNCTION(min) BINARY_FUNCTION(cmptMultiply) BINARY_FUNCTION(cmptDivide) -#undef BINARY_FUNCTION - /* * * * * * * * * * * * * * * * Global operators * * * * * * * * * * * * * */ @@ -143,10 +141,7 @@ inline void opFunc \ UNARY_OPERATOR(-, negate) -#undef UNARY_OPERATOR - - -#define BINARY_OPERATOR_FF(Type1, Type2, op, opFunc) \ +#define BINARY_OPERATOR(Type1, Type2, op, opFunc) \ \ template \ inline void opFunc \ @@ -157,25 +152,11 @@ inline void opFunc \ ) \ {} -#define BINARY_OPERATOR_R(Type1, Type2, op, opFunc) \ - BINARY_OPERATOR_FF(Type1, Type2, op, opFunc) +BINARY_OPERATOR(scalar, Type, *, multiply) +BINARY_OPERATOR(Type, scalar, *, multiply) +BINARY_OPERATOR(Type, scalar, /, divide) -BINARY_OPERATOR_R(Type, Type, +, add) -BINARY_OPERATOR_R(Type, Type, -, subtract) -BINARY_OPERATOR_FF(scalar, Type, *, multiply) -BINARY_OPERATOR_FF(Type, scalar, /, divide) - -#undef BINARY_OPERATOR_R -#undef BINARY_OPERATOR_FF -#undef BINARY_OPERATOR_FTR -#undef BINARY_OPERATOR_TF -#undef BINARY_OPERATOR_TTR -#undef BINARY_OPERATOR_FT -#undef BINARY_OPERATOR_TRF -#undef BINARY_OPERATOR_TRT - - -#define BINARY_TYPE_OPERATOR_TF(TYPE, op, opFunc) \ +#define BINARY_TYPE_OPERATOR_SF(TYPE, op, opFunc) \ \ template \ inline void opFunc \ @@ -187,7 +168,7 @@ inline void opFunc \ {} -#define BINARY_TYPE_OPERATOR_FT(TYPE, op, opFunc) \ +#define BINARY_TYPE_OPERATOR_FS(TYPE, op, opFunc) \ \ template \ inline void opFunc \ @@ -199,19 +180,9 @@ inline void opFunc \ {} -#define BINARY_TYPE_OPERATOR(TYPE, op, opFunc) \ - BINARY_TYPE_OPERATOR_TF(TYPE, op, opFunc) \ - BINARY_TYPE_OPERATOR_FT(TYPE, op, opFunc) - -BINARY_TYPE_OPERATOR(Type, +, add) -BINARY_TYPE_OPERATOR(Type, -, subtract) - -BINARY_TYPE_OPERATOR(scalar, *, multiply) -BINARY_TYPE_OPERATOR_FT(scalar, /, divide) - -#undef BINARY_TYPE_OPERATOR -#undef BINARY_TYPE_OPERATOR_TF -#undef BINARY_TYPE_OPERATOR_FT +BINARY_TYPE_OPERATOR_SF(scalar, *, multiply) +BINARY_TYPE_OPERATOR_FS(scalar, *, multiply) +BINARY_TYPE_OPERATOR_FS(scalar, /, divide) #define PRODUCT_OPERATOR(product, op, opFunc) \ @@ -262,6 +233,9 @@ inline void opFunc \ ) \ {} +PRODUCT_OPERATOR(typeOfSum, +, add) +PRODUCT_OPERATOR(typeOfSum, -, subtract) + PRODUCT_OPERATOR(outerProduct, *, outer) PRODUCT_OPERATOR(crossProduct, ^, cross) PRODUCT_OPERATOR(innerProduct, &, dot) @@ -367,3 +341,7 @@ inline void eigenVectors } // End namespace Foam // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "undefFieldFunctionsM.H" + +// ************************************************************************* // diff --git a/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/tabulated6DoFMotion/tabulated6DoFMotion.C b/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/tabulated6DoFMotion/tabulated6DoFMotion.C index a6924cc4a6..adb53704e9 100644 --- a/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/tabulated6DoFMotion/tabulated6DoFMotion.C +++ b/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/tabulated6DoFMotion/tabulated6DoFMotion.C @@ -127,7 +127,10 @@ bool Foam::solidBodyMotionFunctions::tabulated6DoFMotion::read // If the timeDataFileName has changed read the file - fileName newTimeDataFileName(SBMFCoeffs_.lookup("timeDataFileName")); + fileName newTimeDataFileName + ( + fileName(SBMFCoeffs_.lookup("timeDataFileName")).expand() + ); if (newTimeDataFileName != timeDataFileName_) { diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files index 126669d56a..9c280f9288 100644 --- a/src/finiteVolume/Make/files +++ b/src/finiteVolume/Make/files @@ -158,6 +158,7 @@ $(derivedFvPatchFields)/waveTransmissive/waveTransmissiveFvPatchFields.C $(derivedFvPatchFields)/uniformDensityHydrostaticPressure/uniformDensityHydrostaticPressureFvPatchScalarField.C $(derivedFvPatchFields)/swirlFlowRateInletVelocity/swirlFlowRateInletVelocityFvPatchVectorField.C $(derivedFvPatchFields)/cylindricalInletVelocity/cylindricalInletVelocityFvPatchVectorField.C +$(derivedFvPatchFields)/outletMappedUniformInlet/outletMappedUniformInletFvPatchFields.C fvsPatchFields = fields/fvsPatchFields diff --git a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.C new file mode 100644 index 0000000000..c3f2d89451 --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.C @@ -0,0 +1,190 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ 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 . + +\*---------------------------------------------------------------------------*/ + +#include "outletMappedUniformInletFvPatchField.H" +#include "volFields.H" +#include "surfaceFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +outletMappedUniformInletFvPatchField:: +outletMappedUniformInletFvPatchField +( + const fvPatch& p, + const DimensionedField& iF +) +: + fixedValueFvPatchField(p, iF), + outletPatchName_(), + phiName_("phi") +{} + + +template +outletMappedUniformInletFvPatchField:: +outletMappedUniformInletFvPatchField +( + const outletMappedUniformInletFvPatchField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + fixedValueFvPatchField(ptf, p, iF, mapper), + outletPatchName_(ptf.outletPatchName_), + phiName_(ptf.phiName_) +{} + + +template +outletMappedUniformInletFvPatchField:: +outletMappedUniformInletFvPatchField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + fixedValueFvPatchField(p, iF, dict), + outletPatchName_(dict.lookup("outletPatchName")), + phiName_(dict.lookupOrDefault("phi", "phi")) +{} + + +template +outletMappedUniformInletFvPatchField:: +outletMappedUniformInletFvPatchField +( + const outletMappedUniformInletFvPatchField& ptf +) +: + fixedValueFvPatchField(ptf), + outletPatchName_(ptf.outletPatchName_), + phiName_(ptf.phiName_) +{} + + + +template +outletMappedUniformInletFvPatchField:: +outletMappedUniformInletFvPatchField +( + const outletMappedUniformInletFvPatchField& ptf, + const DimensionedField& iF +) +: + fixedValueFvPatchField(ptf, iF), + outletPatchName_(ptf.outletPatchName_), + phiName_(ptf.phiName_) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +void outletMappedUniformInletFvPatchField::updateCoeffs() +{ + if (this->updated()) + { + return; + } + + const GeometricField& f + ( + dynamic_cast&> + ( + this->dimensionedInternalField() + ) + ); + + const fvPatch& p = this->patch(); + label outletPatchID = + p.patch().boundaryMesh().findPatchID(outletPatchName_); + + if (outletPatchID < 0) + { + FatalErrorIn + ( + "void outletMappedUniformInletFvPatchField::updateCoeffs()" + ) << "Unable to find outlet patch " << outletPatchName_ + << abort(FatalError); + } + + const fvPatch& outletPatch = p.boundaryMesh()[outletPatchID]; + + const fvPatchField& outletPatchField = + f.boundaryField()[outletPatchID]; + + const surfaceScalarField& phi = + this->db().objectRegistry::lookupObject(phiName_); + const scalarField& outletPatchPhi = phi.boundaryField()[outletPatchID]; + scalar sumOutletPatchPhi = gSum(outletPatchPhi); + + if (sumOutletPatchPhi > SMALL) + { + Type averageOutletField = + gSum(outletPatchPhi*outletPatchField) + /sumOutletPatchPhi; + + this->operator==(averageOutletField); + } + else + { + Type averageOutletField = + gSum(outletPatch.magSf()*outletPatchField) + /gSum(outletPatch.magSf()); + + this->operator==(averageOutletField); + } + + fixedValueFvPatchField::updateCoeffs(); +} + + +template +void outletMappedUniformInletFvPatchField::write(Ostream& os) const +{ + fvPatchField::write(os); + os.writeKeyword("outletPatchName") + << outletPatchName_ << token::END_STATEMENT << nl; + if (phiName_ != "phi") + { + os.writeKeyword("phi") << phiName_ << token::END_STATEMENT << nl; + } + this->writeEntry("value", os); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.H b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.H new file mode 100644 index 0000000000..1dee089380 --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.H @@ -0,0 +1,169 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ 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 . + +Class + Foam::outletMappedUniformInletFvPatchField + +Description + Averages the field over the "outlet" patch specified by name + "outletPatchName" and applies this as the uniform value of the field + over this patch. + +SourceFiles + outletMappedUniformInletFvPatchField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef outletMappedUniformInletFvPatchField_H +#define outletMappedUniformInletFvPatchField_H + +#include "fixedValueFvPatchFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class outletMappedUniformInletFvPatch Declaration +\*---------------------------------------------------------------------------*/ + +template +class outletMappedUniformInletFvPatchField +: + public fixedValueFvPatchField +{ + // Private data + + //- Name of the outlet patch to be mapped + word outletPatchName_; + + //- Name of the flux transporting the field + word phiName_; + + +public: + + //- Runtime type information + TypeName("outletMappedUniformInlet"); + + + // Constructors + + //- Construct from patch and internal field + outletMappedUniformInletFvPatchField + ( + const fvPatch&, + const DimensionedField& + ); + + //- Construct from patch, internal field and dictionary + outletMappedUniformInletFvPatchField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given outletMappedUniformInletFvPatchField + // onto a new patch + outletMappedUniformInletFvPatchField + ( + const outletMappedUniformInletFvPatchField&, + const fvPatch&, + const DimensionedField&, + const fvPatchFieldMapper& + ); + + //- Construct as copy + outletMappedUniformInletFvPatchField + ( + const outletMappedUniformInletFvPatchField& + ); + + //- Construct and return a clone + virtual tmp > clone() const + { + return tmp > + ( + new outletMappedUniformInletFvPatchField(*this) + ); + } + + //- Construct as copy setting internal field reference + outletMappedUniformInletFvPatchField + ( + const outletMappedUniformInletFvPatchField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp > clone + ( + const DimensionedField& iF + ) const + { + return tmp > + ( + new outletMappedUniformInletFvPatchField(*this, iF) + ); + } + + + // Member functions + + // Access + + //- Name of the outlet patch to be mapped + const word& outletPatchName() const + { + return outletPatchName_; + } + + + // Evaluation functions + + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + + //- Write + virtual void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "outletMappedUniformInletFvPatchField.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFields.C b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFields.C new file mode 100644 index 0000000000..a932dca978 --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFields.C @@ -0,0 +1,43 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ 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 . + +\*---------------------------------------------------------------------------*/ + +#include "outletMappedUniformInletFvPatchFields.H" +#include "addToRunTimeSelectionTable.H" +#include "volFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +makePatchFields(outletMappedUniformInlet); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFields.H b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFields.H new file mode 100644 index 0000000000..17e96617df --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFields.H @@ -0,0 +1,49 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ 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 . + +\*---------------------------------------------------------------------------*/ + +#ifndef outletMappedUniformInletFvPatchFields_H +#define outletMappedUniformInletFvPatchFields_H + +#include "outletMappedUniformInletFvPatchField.H" +#include "fieldTypes.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +makePatchTypeFieldTypedefs(outletMappedUniformInlet) + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFieldsFwd.H b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFieldsFwd.H new file mode 100644 index 0000000000..76033b9923 --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchFieldsFwd.H @@ -0,0 +1,50 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ 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 . + +\*---------------------------------------------------------------------------*/ + +#ifndef outletMappedUniformInletFvPatchFieldsFwd_H +#define outletMappedUniformInletFvPatchFieldsFwd_H + +#include "fieldTypes.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +template class outletMappedUniformInletFvPatchField; + +makePatchTypeFieldTypedefs(outletMappedUniformInlet) + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C index ea47c95736..ed1e524bf6 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C @@ -423,23 +423,21 @@ void Foam::sixDoFRigidBodyMotion::updateForce scalar deltaT ) { - vector a = vector::zero; + vector fGlobal = vector::zero; - vector tau = vector::zero; + vector tauGlobal = vector::zero; if (Pstream::master()) { + fGlobal = sum(forces); + forAll(positions, i) { - const vector& f = forces[i]; - - a += f/mass_; - - tau += Q().T() & ((positions[i] - centreOfMass()) ^ f); + tauGlobal += (positions[i] - centreOfMass()) ^ forces[i]; } } - updateForce(a, tau, deltaT); + updateForce(fGlobal, tauGlobal, deltaT); }