From e5cd8fff64705628ced1701380b0fff2819d7179 Mon Sep 17 00:00:00 2001 From: Henry Date: Fri, 4 May 2012 16:20:11 +0100 Subject: [PATCH] compressible solvers: renamed Kp -> Ekp (kinetic + pressure energy) --- .../solvers/compressible/rhoSimpleFoam/eEqn.H | 18 ++++++++++++++++++ .../rhoPorousMRFSimpleFoam/eEqn.H | 19 +++++++++++++++++++ .../angledDuctExplicit/system/fvSchemes | 6 +++--- .../angledDuctExplicit/system/fvSolution | 4 ++-- .../constant/thermophysicalProperties | 4 ++-- .../angledDuctImplicit/system/fvSchemes | 6 +++--- .../angledDuctImplicit/system/fvSolution | 4 ++-- .../constant/thermophysicalProperties | 4 ++-- .../squareBend/system/fvSchemes | 5 +++-- .../squareBend/system/fvSolution | 4 ++-- 10 files changed, 56 insertions(+), 18 deletions(-) create mode 100644 applications/solvers/compressible/rhoSimpleFoam/eEqn.H create mode 100644 applications/solvers/compressible/rhoSimpleFoam/rhoPorousMRFSimpleFoam/eEqn.H diff --git a/applications/solvers/compressible/rhoSimpleFoam/eEqn.H b/applications/solvers/compressible/rhoSimpleFoam/eEqn.H new file mode 100644 index 0000000000..a1ea771573 --- /dev/null +++ b/applications/solvers/compressible/rhoSimpleFoam/eEqn.H @@ -0,0 +1,18 @@ +{ + // Kinetic + pressure energy + volScalarField Ekp("Ekp", 0.5*magSqr(U) + p/rho); + + fvScalarMatrix eEqn + ( + fvm::div(phi, e) + - fvm::Sp(fvc::div(phi), e) + - fvm::laplacian(turbulence->alphaEff(), e) + == + fvc::div(phi)*Ekp - fvc::div(phi, Ekp) + ); + + eEqn.relax(); + eEqn.solve(); + + thermo.correct(); +} diff --git a/applications/solvers/compressible/rhoSimpleFoam/rhoPorousMRFSimpleFoam/eEqn.H b/applications/solvers/compressible/rhoSimpleFoam/rhoPorousMRFSimpleFoam/eEqn.H new file mode 100644 index 0000000000..4791062d35 --- /dev/null +++ b/applications/solvers/compressible/rhoSimpleFoam/rhoPorousMRFSimpleFoam/eEqn.H @@ -0,0 +1,19 @@ +{ + volScalarField Ekp("Ekp", 0.5*magSqr(U) + p/rho); + + fvScalarMatrix eEqn + ( + fvm::div(phi, e) + - fvm::Sp(fvc::div(phi), e) + - fvm::laplacian(turbulence->alphaEff(), e) + == + fvc::div(phi)*Ekp - fvc::div(phi, Ekp) + ); + + //pZones.addEnergySource(thermo, rho, eEqn); + + eEqn.relax(); + eEqn.solve(); + + thermo.correct(); +} diff --git a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSchemes b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSchemes index 711e4e3233..0a14cc31a5 100644 --- a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSchemes +++ b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSchemes @@ -29,16 +29,16 @@ divSchemes { div(phi,U) Gauss upwind; div((muEff*dev2(T(grad(U))))) Gauss linear; - div(phi,h) Gauss upwind; + div(phi,e) Gauss upwind; div(phi,epsilon) Gauss upwind; div(phi,k) Gauss upwind; - div(phi,K) Gauss upwind; + div(phi,Ekp) Gauss upwind; } laplacianSchemes { laplacian(muEff,U) Gauss linear corrected; - laplacian(alphaEff,h) Gauss linear corrected; + laplacian(alphaEff,e) Gauss linear corrected; laplacian((rho*rAU),p) Gauss linear corrected; laplacian(DepsilonEff,epsilon) Gauss linear corrected; laplacian(DkEff,k) Gauss linear corrected; diff --git a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSolution b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSolution index 317c954ef6..a1e3109222 100644 --- a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSolution +++ b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctExplicit/system/fvSolution @@ -38,7 +38,7 @@ solvers relTol 0.1; } - h + e { solver PBiCG; preconditioner DILU; @@ -82,7 +82,7 @@ relaxationFactors { U 0.7; "(k|epsilon)" 0.7; - h 0.5; + e 0.5; } } diff --git a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/constant/thermophysicalProperties b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/constant/thermophysicalProperties index d6d597d433..cbe2156297 100644 --- a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/constant/thermophysicalProperties +++ b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/constant/thermophysicalProperties @@ -15,7 +15,7 @@ FoamFile } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -thermoType hPsiThermo>>>>; +thermoType ePsiThermo>>>>; mixture { @@ -26,7 +26,7 @@ mixture } thermodynamics { - Cp 1007; + Cv 719.3; Hf 0; } transport diff --git a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSchemes b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSchemes index 711e4e3233..0a14cc31a5 100644 --- a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSchemes +++ b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSchemes @@ -29,16 +29,16 @@ divSchemes { div(phi,U) Gauss upwind; div((muEff*dev2(T(grad(U))))) Gauss linear; - div(phi,h) Gauss upwind; + div(phi,e) Gauss upwind; div(phi,epsilon) Gauss upwind; div(phi,k) Gauss upwind; - div(phi,K) Gauss upwind; + div(phi,Ekp) Gauss upwind; } laplacianSchemes { laplacian(muEff,U) Gauss linear corrected; - laplacian(alphaEff,h) Gauss linear corrected; + laplacian(alphaEff,e) Gauss linear corrected; laplacian((rho*rAU),p) Gauss linear corrected; laplacian(DepsilonEff,epsilon) Gauss linear corrected; laplacian(DkEff,k) Gauss linear corrected; diff --git a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSolution b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSolution index 6f93606315..65a9606e7f 100644 --- a/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSolution +++ b/tutorials/compressible/rhoPorousMRFSimpleFoam/angledDuctImplicit/system/fvSolution @@ -29,7 +29,7 @@ solvers mergeLevels 1; } - h + e { solver PBiCG; preconditioner DILU; @@ -76,7 +76,7 @@ relaxationFactors { U 0.7; "(k|epsilon)" 0.9; - h 0.9; + e 0.9; } } diff --git a/tutorials/compressible/rhoSimplecFoam/squareBend/constant/thermophysicalProperties b/tutorials/compressible/rhoSimplecFoam/squareBend/constant/thermophysicalProperties index d6d597d433..cbe2156297 100644 --- a/tutorials/compressible/rhoSimplecFoam/squareBend/constant/thermophysicalProperties +++ b/tutorials/compressible/rhoSimplecFoam/squareBend/constant/thermophysicalProperties @@ -15,7 +15,7 @@ FoamFile } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -thermoType hPsiThermo>>>>; +thermoType ePsiThermo>>>>; mixture { @@ -26,7 +26,7 @@ mixture } thermodynamics { - Cp 1007; + Cv 719.3; Hf 0; } transport diff --git a/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSchemes b/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSchemes index c9c95b180c..e1e47211d7 100644 --- a/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSchemes +++ b/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSchemes @@ -31,12 +31,13 @@ divSchemes div(phi,U) Gauss upwind; div((muEff*dev2(T(grad(U))))) Gauss linear; - div(phi,h) Gauss upwind; + div(phi,e) Gauss upwind; div(phi,epsilon) Gauss upwind; div(phi,k) Gauss upwind; div(phid,p) Gauss upwind; - div(phi,K) Gauss upwind; + div(phi,Ekp) Gauss upwind; + div((phi|interpolate(rho)),p) Gauss upwind; } laplacianSchemes diff --git a/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSolution b/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSolution index 09afca5761..812766041a 100644 --- a/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSolution +++ b/tutorials/compressible/rhoSimplecFoam/squareBend/system/fvSolution @@ -32,7 +32,7 @@ solvers mergeLevels 1; } - "(U|h|k|epsilon)" + "(U|e|k|epsilon)" { solver GAMG; tolerance 1e-08; @@ -67,7 +67,7 @@ relaxationFactors { p 1; U 0.9; - h 0.9; + e 0.9; k 0.9; epsilon 0.9; }