From 7b4002fec2abb02fdb1ad5d48381843a016635ef Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Sat, 25 Mar 2023 15:36:09 +0000 Subject: [PATCH] solvers::*:moveMesh: rationalised the application of correctPhi Flux correction is now applied if either the topology changed or the mesh is moving and correctPhi is true. This strategy allows moving mesh cases without topology change to be run without any change to the fluxes which is appropriate for solid-body motion of the entire domain or a rotating sub-domain with NCC. --- .../solvers/modules/VoFSolver/moveMesh.C | 53 ++++++++-------- .../modules/incompressibleFluid/moveMesh.C | 23 ++++--- .../modules/isothermalFluid/moveMesh.C | 23 ++++--- .../multiphaseEuler/moveMesh.C | 28 +++++---- .../phaseSystems/phaseSystem/phaseSystem.C | 60 +++++++++---------- .../phaseSystems/phaseSystem/phaseSystem.H | 1 - 6 files changed, 89 insertions(+), 99 deletions(-) diff --git a/applications/solvers/modules/VoFSolver/moveMesh.C b/applications/solvers/modules/VoFSolver/moveMesh.C index f3e9248781..a61004238c 100644 --- a/applications/solvers/modules/VoFSolver/moveMesh.C +++ b/applications/solvers/modules/VoFSolver/moveMesh.C @@ -35,7 +35,7 @@ void Foam::solvers::VoFSolver::moveMesh() { if ( - correctPhi + (correctPhi || mesh.topoChanged()) && divergent() && !divU.valid() ) @@ -65,34 +65,31 @@ void Foam::solvers::VoFSolver::moveMesh() correctUphiBCs(U_, phi, true); - if (correctPhi) + if (divU.valid()) { - if (divU.valid()) - { - CorrectPhi - ( - phi, - U, - p_rgh, - surfaceScalarField("rAUf", fvc::interpolate(rAU())), - divU(), - pressureReference(), - pimple - ); - } - else - { - CorrectPhi - ( - phi, - U, - p_rgh, - surfaceScalarField("rAUf", fvc::interpolate(rAU())), - geometricZeroField(), - pressureReference(), - pimple - ); - } + CorrectPhi + ( + phi, + U, + p_rgh, + surfaceScalarField("rAUf", fvc::interpolate(rAU())), + divU(), + pressureReference(), + pimple + ); + } + else + { + CorrectPhi + ( + phi, + U, + p_rgh, + surfaceScalarField("rAUf", fvc::interpolate(rAU())), + geometricZeroField(), + pressureReference(), + pimple + ); } // Make the fluxes relative to the mesh motion diff --git a/applications/solvers/modules/incompressibleFluid/moveMesh.C b/applications/solvers/modules/incompressibleFluid/moveMesh.C index 0659378aa3..6d34cb773d 100644 --- a/applications/solvers/modules/incompressibleFluid/moveMesh.C +++ b/applications/solvers/modules/incompressibleFluid/moveMesh.C @@ -48,19 +48,16 @@ void Foam::solvers::incompressibleFluid::moveMesh() correctUphiBCs(U, phi, true); - if (correctPhi) - { - CorrectPhi - ( - phi, - U, - p, - dimensionedScalar("rAUf", dimTime, 1), - geometricZeroField(), - pressureReference, - pimple - ); - } + 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 af63d3d129..b5ebf5c291 100644 --- a/applications/solvers/modules/isothermalFluid/moveMesh.C +++ b/applications/solvers/modules/isothermalFluid/moveMesh.C @@ -57,19 +57,16 @@ void Foam::solvers::isothermalFluid::moveMesh() correctUphiBCs(rho, U, phi, true); - if (correctPhi) - { - CorrectPhi - ( - phi, - buoyancy.valid() ? p_rgh : p, - rho, - thermo.psi(), - dimensionedScalar("rAUf", dimTime, 1), - divrhoU(), - pimple - ); - } + 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 9379062b4e..a843cc4de9 100644 --- a/applications/solvers/modules/multiphaseEuler/multiphaseEuler/moveMesh.C +++ b/applications/solvers/modules/multiphaseEuler/multiphaseEuler/moveMesh.C @@ -37,7 +37,11 @@ void Foam::solvers::multiphaseEuler::moveMesh() && (pimple.firstIter() || pimple.moveMeshOuterCorrectors()) ) { - if (correctPhi && !divU.valid()) + if + ( + (correctPhi || mesh.topoChanged()) + && !divU.valid() + ) { // Construct and register divU for mapping divU = new volScalarField @@ -53,20 +57,22 @@ void Foam::solvers::multiphaseEuler::moveMesh() // Move the mesh mesh_.move(); - if (mesh.changing() || mesh.topoChanged()) + if (mesh.changing()) { buoyancy.moveMesh(); - fluid.meshUpdate(); + if (correctPhi || mesh.topoChanged()) + { + fluid.meshUpdate(); - fluid.correctPhi - ( - p_rgh, - divU, - correctPhi, - pressureReference, - pimple - ); + fluid.correctPhi + ( + p_rgh, + divU, + pressureReference, + pimple + ); + } meshCourantNo(); } diff --git a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C index 39e2859fc2..de3c973dee 100644 --- a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C +++ b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.C @@ -727,7 +727,6 @@ void Foam::phaseSystem::correctPhi ( const volScalarField& p_rgh, const tmp& divU, - const bool correctPhi, const pressureReference& pressureReference, nonOrthogonalSolutionControl& pimple ) @@ -760,42 +759,37 @@ void Foam::phaseSystem::correctPhi // Correct fixed-flux BCs to be consistent with the velocity BCs correctBoundaryFlux(); + phi_ = Zero; + PtrList alphafs(phaseModels_.size()); + forAll(movingPhases(), movingPhasei) { - phi_ = Zero; - PtrList alphafs(phaseModels_.size()); - forAll(movingPhases(), movingPhasei) - { - phaseModel& phase = movingPhases()[movingPhasei]; - const label phasei = phase.index(); - const volScalarField& alpha = phase; + phaseModel& phase = movingPhases()[movingPhasei]; + const label phasei = phase.index(); + const volScalarField& alpha = phase; - alphafs.set(phasei, fvc::interpolate(alpha).ptr()); + alphafs.set(phasei, fvc::interpolate(alpha).ptr()); - // Calculate absolute flux - // from the mapped surface velocity - phi_ += alphafs[phasei]*(mesh_.Sf() & phase.UfRef()); - } - - 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()); - - setMixturePhi(alphafs, phi_); + // Calculate absolute flux + // from the mapped surface velocity + phi_ += alphafs[phasei]*(mesh_.Sf() & phase.UfRef()); } + + 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()); + + setMixturePhi(alphafs, phi_); } diff --git a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H index 035c52050d..09a00702cb 100644 --- a/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H +++ b/applications/solvers/modules/multiphaseEuler/phaseSystems/phaseSystem/phaseSystem.H @@ -657,7 +657,6 @@ public: ( const volScalarField& p_rgh, const tmp& divU, - const bool correctPhi, const pressureReference& pressureReference, nonOrthogonalSolutionControl& pimple );