From 084097bab9299f2417652dd744d0c5cd4d22d551 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Sun, 20 Nov 2022 20:22:59 +0000 Subject: [PATCH] solvers/modules: Updated correctPhi logic so that it is not required following topology change Now fluxes are updated from the mapped fields following mesh topology change with or without implicit continuity correction enabled by the optional correctPhi switch. --- .../modules/compressibleVoF/moveMesh.C | 25 ++++++++++-------- .../modules/incompressibleFluid/moveMesh.C | 25 ++++++++++-------- .../modules/isothermalFluid/moveMesh.C | 25 ++++++++++-------- .../multiphaseEuler/moveMesh.C | 20 +++++++------- .../multiphaseEuler/multiphaseEuler.C | 2 +- .../phaseSystems/phaseSystem/phaseSystem.C | 26 +++++++++++-------- .../phaseSystems/phaseSystem/phaseSystem.H | 1 + 7 files changed, 68 insertions(+), 56 deletions(-) diff --git a/applications/solvers/modules/compressibleVoF/moveMesh.C b/applications/solvers/modules/compressibleVoF/moveMesh.C index 1c22912427..a2655d1df5 100644 --- a/applications/solvers/modules/compressibleVoF/moveMesh.C +++ b/applications/solvers/modules/compressibleVoF/moveMesh.C @@ -45,7 +45,7 @@ bool Foam::solvers::compressibleVoF::moveMesh() // Move the mesh mesh.move(); - if (mesh.changing()) + if (correctPhi || mesh.topoChanged()) { buoyancy.moveMesh(); @@ -59,16 +59,19 @@ bool Foam::solvers::compressibleVoF::moveMesh() correctUphiBCs(U, phi, true); - CorrectPhi - ( - phi, - U, - p_rgh, - surfaceScalarField("rAUf", fvc::interpolate(rAU())), - divU(), - pressureReference, - pimple - ); + if (correctPhi) + { + CorrectPhi + ( + phi, + U, + p_rgh, + surfaceScalarField("rAUf", fvc::interpolate(rAU())), + divU(), + pressureReference, + pimple + ); + } // Make the fluxes relative to the mesh motion fvc::makeRelative(phi, U); diff --git a/applications/solvers/modules/incompressibleFluid/moveMesh.C b/applications/solvers/modules/incompressibleFluid/moveMesh.C index 06761f1349..0357e321e6 100644 --- a/applications/solvers/modules/incompressibleFluid/moveMesh.C +++ b/applications/solvers/modules/incompressibleFluid/moveMesh.C @@ -39,7 +39,7 @@ bool Foam::solvers::incompressibleFluid::moveMesh() { MRF.update(); - if (correctPhi) + if (correctPhi || mesh.topoChanged()) { // Calculate absolute flux // from the mapped surface velocity @@ -47,16 +47,19 @@ bool Foam::solvers::incompressibleFluid::moveMesh() correctUphiBCs(U, phi, true); - CorrectPhi - ( - phi, - U, - p, - dimensionedScalar("rAUf", dimTime, 1), - geometricZeroField(), - pressureReference, - pimple - ); + if (correctPhi) + { + CorrectPhi + ( + phi, + U, + p, + dimensionedScalar("rAUf", dimTime, 1), + geometricZeroField(), + pressureReference, + pimple + ); + } // Make the flux relative to the mesh motion fvc::makeRelative(phi, U); diff --git a/applications/solvers/modules/isothermalFluid/moveMesh.C b/applications/solvers/modules/isothermalFluid/moveMesh.C index d6dc2b9aaa..5e422dc5c5 100644 --- a/applications/solvers/modules/isothermalFluid/moveMesh.C +++ b/applications/solvers/modules/isothermalFluid/moveMesh.C @@ -44,7 +44,7 @@ bool Foam::solvers::isothermalFluid::moveMesh() MRF.update(); - if (correctPhi) + if (correctPhi || mesh.topoChanged()) { // Calculate absolute flux // from the mapped surface velocity @@ -52,16 +52,19 @@ bool Foam::solvers::isothermalFluid::moveMesh() correctUphiBCs(rho, U, phi, true); - CorrectPhi - ( - phi, - buoyancy.valid() ? p_rgh : p, - rho, - thermo.psi(), - dimensionedScalar("rAUf", dimTime, 1), - divrhoU(), - pimple - ); + if (correctPhi) + { + CorrectPhi + ( + phi, + buoyancy.valid() ? p_rgh : p, + rho, + thermo.psi(), + dimensionedScalar("rAUf", dimTime, 1), + divrhoU(), + pimple + ); + } // Make the fluxes relative to the mesh-motion fvc::makeRelative(phi, rho, U); diff --git a/applications/solvers/modules/multiphaseEuler/multiphaseEuler/moveMesh.C b/applications/solvers/modules/multiphaseEuler/multiphaseEuler/moveMesh.C index 4a2c536946..2a79c2ea9f 100644 --- a/applications/solvers/modules/multiphaseEuler/multiphaseEuler/moveMesh.C +++ b/applications/solvers/modules/multiphaseEuler/multiphaseEuler/moveMesh.C @@ -53,22 +53,20 @@ bool Foam::solvers::multiphaseEuler::moveMesh() // Move the mesh mesh.move(); - if (mesh.changing()) + if (mesh.changing() || mesh.topoChanged()) { buoyancy.moveMesh(); fluid.meshUpdate(); - if (correctPhi) - { - fluid.correctPhi - ( - p_rgh, - divU, - pressureReference, - pimple - ); - } + fluid.correctPhi + ( + p_rgh, + divU, + correctPhi, + pressureReference, + pimple + ); meshCourantNo(); diff --git a/applications/solvers/modules/multiphaseEuler/multiphaseEuler/multiphaseEuler.C b/applications/solvers/modules/multiphaseEuler/multiphaseEuler/multiphaseEuler.C index 0203b3f2ba..238385036e 100644 --- a/applications/solvers/modules/multiphaseEuler/multiphaseEuler/multiphaseEuler.C +++ b/applications/solvers/modules/multiphaseEuler/multiphaseEuler/multiphaseEuler.C @@ -207,7 +207,7 @@ void Foam::solvers::multiphaseEuler::preSolve() // Store divU from the previous mesh so that it can be // mapped and used in correctPhi to ensure the corrected phi // has the same divergence - if (correctPhi && mesh.topoChanged()) + if (correctPhi) { // Construct and register divU for mapping divU = new volScalarField diff --git a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C index 52d9c0c20b..81d1298f77 100644 --- a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C +++ b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C @@ -710,6 +710,7 @@ void Foam::phaseSystem::correctPhi ( const volScalarField& p_rgh, const tmp& divU, + const bool correctPhi, const pressureReference& pressureReference, nonOrthogonalSolutionControl& pimple ) @@ -758,17 +759,20 @@ void Foam::phaseSystem::correctPhi phi_ += alphafs[phasei]*(mesh_.Sf() & phase.Uf()); } - CorrectPhi - ( - phi_, - movingPhases()[0].U(), - p_rgh, - // surfaceScalarField("rAUf", fvc::interpolate(rAU())), - dimensionedScalar(dimTime/dimDensity, 1), - divU(), - pressureReference, - pimple - ); + if (correctPhi) + { + CorrectPhi + ( + phi_, + movingPhases()[0].U(), + p_rgh, + // surfaceScalarField("rAUf", fvc::interpolate(rAU())), + dimensionedScalar(dimTime/dimDensity, 1), + divU(), + pressureReference, + pimple + ); + } // Make the flux relative to the mesh motion fvc::makeRelative(phi_, movingPhases()[0].U()); diff --git a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H index d042116f1e..6d4662e746 100644 --- a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H +++ b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H @@ -651,6 +651,7 @@ public: ( const volScalarField& p_rgh, const tmp& divU, + const bool correctPhi, const pressureReference& pressureReference, nonOrthogonalSolutionControl& pimple );