From 8fa6bfcded9cdc3845a6ab6eddbabcec5bf058da Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Thu, 1 Oct 2020 10:44:36 +0100 Subject: [PATCH] compressibleInterFoam: Updated to use the thermo:rho --- .../compressibleAlphaEqnSubCycle.H | 2 +- .../compressibleInterFilmFoam.C | 6 +----- .../compressibleInterFilmFoam/pEqn.H | 3 +-- .../compressibleInterFoam.C | 6 +----- .../compressibleInterFoam/createFieldRefs.H | 9 ++++++++ .../compressibleInterFoam/createFields.H | 21 +------------------ .../multiphase/compressibleInterFoam/pEqn.H | 3 +-- .../laminar/climbingRod/0/p_rgh | 1 + .../laminar/depthCharge2D/system/fvSchemes | 2 +- .../laminar/depthCharge3D/system/fvSchemes | 2 +- .../laminar/sloshingTank2D/system/fvSchemes | 2 +- 11 files changed, 19 insertions(+), 38 deletions(-) create mode 100644 applications/solvers/multiphase/compressibleInterFoam/createFieldRefs.H diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleAlphaEqnSubCycle.H b/applications/solvers/multiphase/compressibleInterFoam/compressibleAlphaEqnSubCycle.H index f7cf0188d9..b48d0ba436 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/compressibleAlphaEqnSubCycle.H +++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleAlphaEqnSubCycle.H @@ -54,7 +54,7 @@ else #include "alphaEqn.H" } -rho == alpha1*rho1 + alpha2*rho2; +// rho == alpha1*rho1 + alpha2*rho2; const surfaceScalarField& alphaPhi1 = talphaPhi1(); surfaceScalarField alphaPhi2("alphaPhi2", phi - alphaPhi1); diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/compressibleInterFilmFoam.C b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/compressibleInterFilmFoam.C index 87f386a471..5103a0e4f6 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/compressibleInterFilmFoam.C +++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/compressibleInterFilmFoam.C @@ -61,13 +61,9 @@ int main(int argc, char *argv[]) #include "createControl.H" #include "createTimeControls.H" #include "createFields.H" + #include "createFieldRefs.H" #include "createSurfaceFilmModel.H" - volScalarField& p = mixture.p(); - volScalarField& T = mixture.T(); - const volScalarField& psi1 = mixture.thermo1().psi(); - const volScalarField& psi2 = mixture.thermo2().psi(); - regionModels::surfaceFilmModel& surfaceFilm = tsurfaceFilm(); if (!LTS) diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/pEqn.H b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/pEqn.H index 7b89c120b9..0535d15d6b 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/pEqn.H +++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFilmFoam/pEqn.H @@ -133,8 +133,7 @@ // Update densities from change in p_rgh mixture.thermo1().correctRho(psi1*(p_rgh - p_rgh_0)); mixture.thermo2().correctRho(psi2*(p_rgh - p_rgh_0)); - - rho = alpha1*rho1 + alpha2*rho2; + mixture.correct(); // Correct p_rgh for consistency with p and the updated densities p_rgh = p - rho*gh; diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFoam.C b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFoam.C index 72e9117e5e..8fcd8de136 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFoam.C +++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterFoam.C @@ -65,15 +65,11 @@ int main(int argc, char *argv[]) #include "initContinuityErrs.H" #include "createDyMControls.H" #include "createFields.H" + #include "createFieldRefs.H" #include "CourantNo.H" #include "setInitialDeltaT.H" #include "createUfIfPresent.H" - volScalarField& p = mixture.p(); - volScalarField& T = mixture.T(); - const volScalarField& psi1 = mixture.thermo1().psi(); - const volScalarField& psi2 = mixture.thermo2().psi(); - // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // Info<< "\nStarting time loop\n" << endl; diff --git a/applications/solvers/multiphase/compressibleInterFoam/createFieldRefs.H b/applications/solvers/multiphase/compressibleInterFoam/createFieldRefs.H new file mode 100644 index 0000000000..7d6360ce52 --- /dev/null +++ b/applications/solvers/multiphase/compressibleInterFoam/createFieldRefs.H @@ -0,0 +1,9 @@ +volScalarField& alpha2(mixture.alpha2()); + +const volScalarField& rho1 = mixture.thermo1().rho(); +const volScalarField& rho2 = mixture.thermo2().rho(); + +volScalarField& p = mixture.p(); +volScalarField& T = mixture.T(); +const volScalarField& psi1 = mixture.thermo1().psi(); +const volScalarField& psi2 = mixture.thermo2().psi(); diff --git a/applications/solvers/multiphase/compressibleInterFoam/createFields.H b/applications/solvers/multiphase/compressibleInterFoam/createFields.H index 5772b179b2..35f6c57ff8 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/createFields.H +++ b/applications/solvers/multiphase/compressibleInterFoam/createFields.H @@ -34,26 +34,7 @@ Info<< "Constructing twoPhaseMixtureThermo\n" << endl; twoPhaseMixtureThermo mixture(U, phi); volScalarField& alpha1(mixture.alpha1()); -volScalarField& alpha2(mixture.alpha2()); - -Info<< "Reading thermophysical properties\n" << endl; - -const volScalarField& rho1 = mixture.thermo1().rho(); -const volScalarField& rho2 = mixture.thermo2().rho(); - -volScalarField rho -( - IOobject - ( - "rho", - runTime.timeName(), - mesh, - IOobject::READ_IF_PRESENT, - IOobject::AUTO_WRITE - ), - alpha1*rho1 + alpha2*rho2 -); - +const volScalarField& rho = mixture.rho(); dimensionedScalar pMin ( diff --git a/applications/solvers/multiphase/compressibleInterFoam/pEqn.H b/applications/solvers/multiphase/compressibleInterFoam/pEqn.H index df0c29fc72..a9535a3e21 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/pEqn.H +++ b/applications/solvers/multiphase/compressibleInterFoam/pEqn.H @@ -150,8 +150,7 @@ // Update densities from change in p_rgh mixture.thermo1().correctRho(psi1*(p_rgh - p_rgh_0)); mixture.thermo2().correctRho(psi2*(p_rgh - p_rgh_0)); - - rho = alpha1*rho1 + alpha2*rho2; + mixture.correct(); // Correct p_rgh for consistency with p and the updated densities p_rgh = p - rho*gh; diff --git a/tutorials/multiphase/compressibleInterFoam/laminar/climbingRod/0/p_rgh b/tutorials/multiphase/compressibleInterFoam/laminar/climbingRod/0/p_rgh index 5f7d8b2839..a7d26832d8 100644 --- a/tutorials/multiphase/compressibleInterFoam/laminar/climbingRod/0/p_rgh +++ b/tutorials/multiphase/compressibleInterFoam/laminar/climbingRod/0/p_rgh @@ -30,6 +30,7 @@ boundaryField { type totalPressure; p0 $internalField; + rho thermo:rho; } #includeEtc "caseDicts/setConstraintTypes" diff --git a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSchemes b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSchemes index 8a496ce3d0..a815750fa7 100644 --- a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSchemes +++ b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSchemes @@ -35,7 +35,7 @@ divSchemes div(phi,p) Gauss upwind; div(phi,k) Gauss upwind; - div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; + div(((thermo:rho*nuEff)*dev2(T(grad(U))))) Gauss linear; } laplacianSchemes diff --git a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSchemes b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSchemes index 8a496ce3d0..a815750fa7 100644 --- a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSchemes +++ b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSchemes @@ -35,7 +35,7 @@ divSchemes div(phi,p) Gauss upwind; div(phi,k) Gauss upwind; - div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; + div(((thermo:rho*nuEff)*dev2(T(grad(U))))) Gauss linear; } laplacianSchemes diff --git a/tutorials/multiphase/compressibleInterFoam/laminar/sloshingTank2D/system/fvSchemes b/tutorials/multiphase/compressibleInterFoam/laminar/sloshingTank2D/system/fvSchemes index 4dfec4d3ed..7b6058c3ae 100644 --- a/tutorials/multiphase/compressibleInterFoam/laminar/sloshingTank2D/system/fvSchemes +++ b/tutorials/multiphase/compressibleInterFoam/laminar/sloshingTank2D/system/fvSchemes @@ -34,7 +34,7 @@ divSchemes div(rhoPhi,K) Gauss linear; div((phi+meshPhi),p) Gauss linear; - div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; + div(((thermo:rho*nuEff)*dev2(T(grad(U))))) Gauss linear; } laplacianSchemes