From eb3e43d1afa9971119cd3a0cb2d572d4dda2e897 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Tue, 17 Sep 2013 16:04:22 +0100 Subject: [PATCH] Moved solvers out of tutorials --- .../incompressible/MRFSimpleFoam/MRFSimpleFoam.C | 129 ++++++++++++++++++++ .../incompressible/MRFSimpleFoam/Make/files | 3 + .../incompressible/MRFSimpleFoam/Make/options | 10 ++ .../incompressible/MRFSimpleFoam/createFields.H | 46 ++++++++ .../incompressible/simpleSRFFoam/Make/files | 3 + .../incompressible/simpleSRFFoam/Make/options | 11 ++ .../incompressible/simpleSRFFoam/createFields.H | 74 ++++++++++++ .../incompressible/simpleSRFFoam/simpleSRFFoam.C | 131 +++++++++++++++++++++ .../solvers/multiphase/MRFInterFoam/MRFInterFoam.C | 109 +++++++++++++++++ .../solvers/multiphase/MRFInterFoam/Make/files | 3 + .../solvers/multiphase/MRFInterFoam/Make/options | 15 +++ .../solvers/multiphase/MRFInterFoam/UEqn.H | 35 ++++++ .../multiphase/MRFInterFoam/createMRFZones.H | 2 + .../solvers/multiphase/MRFInterFoam/pEqn.H | 49 ++++++++ 14 files changed, 620 insertions(+) create mode 100644 applications/solvers/incompressible/MRFSimpleFoam/MRFSimpleFoam.C create mode 100644 applications/solvers/incompressible/MRFSimpleFoam/Make/files create mode 100644 applications/solvers/incompressible/MRFSimpleFoam/Make/options create mode 100644 applications/solvers/incompressible/MRFSimpleFoam/createFields.H create mode 100644 applications/solvers/incompressible/simpleSRFFoam/Make/files create mode 100644 applications/solvers/incompressible/simpleSRFFoam/Make/options create mode 100644 applications/solvers/incompressible/simpleSRFFoam/createFields.H create mode 100644 applications/solvers/incompressible/simpleSRFFoam/simpleSRFFoam.C create mode 100644 applications/solvers/multiphase/MRFInterFoam/MRFInterFoam.C create mode 100644 applications/solvers/multiphase/MRFInterFoam/Make/files create mode 100644 applications/solvers/multiphase/MRFInterFoam/Make/options create mode 100644 applications/solvers/multiphase/MRFInterFoam/UEqn.H create mode 100644 applications/solvers/multiphase/MRFInterFoam/createMRFZones.H create mode 100644 applications/solvers/multiphase/MRFInterFoam/pEqn.H diff --git a/applications/solvers/incompressible/MRFSimpleFoam/MRFSimpleFoam.C b/applications/solvers/incompressible/MRFSimpleFoam/MRFSimpleFoam.C new file mode 100644 index 000000000..4e741c4a8 --- /dev/null +++ b/applications/solvers/incompressible/MRFSimpleFoam/MRFSimpleFoam.C @@ -0,0 +1,129 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright held by original author + \\/ 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 2 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, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Application + MRFSimpleFoam + +Description + Steady-state solver for incompressible, turbulent flow of non-Newtonian + fluids with MRF regions. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "incompressible/singlePhaseTransportModel/singlePhaseTransportModel.H" +#include "incompressible/RAS/RASModel/RASModel.H" +#include "MRFZones.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + +# include "setRootCase.H" + +# include "createTime.H" +# include "createMesh.H" +# include "createFields.H" +# include "initContinuityErrs.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + for (runTime++; !runTime.end(); runTime++) + { + Info<< "Time = " << runTime.timeName() << nl << endl; + +# include "readSIMPLEControls.H" + + p.storePrevIter(); + + // Pressure-velocity SIMPLE corrector + { + // Momentum predictor + + tmp UEqn + ( + fvm::div(phi, U) + + turbulence->divDevReff(U) + ); + mrfZones.addCoriolis(UEqn()); + + UEqn().relax(); + + solve(UEqn() == -fvc::grad(p)); + + p.boundaryField().updateCoeffs(); + volScalarField rAU = 1.0/UEqn().A(); + U = rAU*UEqn().H(); + UEqn.clear(); + + phi = fvc::interpolate(U, "interpolate(HbyA)") & mesh.Sf(); + mrfZones.relativeFlux(phi); + adjustPhi(phi, U, p); + + // Non-orthogonal pressure corrector loop + for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + { + fvScalarMatrix pEqn + ( + fvm::laplacian(rAU, p) == fvc::div(phi) + ); + + pEqn.setReference(pRefCell, pRefValue); + pEqn.solve(); + + if (nonOrth == nNonOrthCorr) + { + phi -= pEqn.flux(); + } + } + +# include "continuityErrs.H" + + // Explicitly relax pressure for momentum corrector + p.relax(); + + // Momentum corrector + U -= rAU*fvc::grad(p); + U.correctBoundaryConditions(); + } + + turbulence->correct(); + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + } + + Info<< "End\n" << endl; + + return(0); +} + + +// ************************************************************************* // diff --git a/applications/solvers/incompressible/MRFSimpleFoam/Make/files b/applications/solvers/incompressible/MRFSimpleFoam/Make/files new file mode 100644 index 000000000..56c194272 --- /dev/null +++ b/applications/solvers/incompressible/MRFSimpleFoam/Make/files @@ -0,0 +1,3 @@ +MRFSimpleFoam.C + +EXE = $(FOAM_APPBIN)/MRFSimpleFoam diff --git a/applications/solvers/incompressible/MRFSimpleFoam/Make/options b/applications/solvers/incompressible/MRFSimpleFoam/Make/options new file mode 100644 index 000000000..44cec8365 --- /dev/null +++ b/applications/solvers/incompressible/MRFSimpleFoam/Make/options @@ -0,0 +1,10 @@ +EXE_INC = \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/turbulenceModels \ + -I$(LIB_SRC)/transportModels + +EXE_LIBS = \ + -lincompressibleRASModels \ + -lincompressibleTransportModels \ + -lfiniteVolume \ + -llduSolvers diff --git a/applications/solvers/incompressible/MRFSimpleFoam/createFields.H b/applications/solvers/incompressible/MRFSimpleFoam/createFields.H new file mode 100644 index 000000000..43242a6e3 --- /dev/null +++ b/applications/solvers/incompressible/MRFSimpleFoam/createFields.H @@ -0,0 +1,46 @@ + Info << "Reading field p\n" << endl; + volScalarField p + ( + IOobject + ( + "p", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + Info << "Reading field U\n" << endl; + volVectorField U + ( + IOobject + ( + "U", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + +# include "createPhi.H" + + + label pRefCell = 0; + scalar pRefValue = 0.0; + setRefCell(p, mesh.solutionDict().subDict("SIMPLE"), pRefCell, pRefValue); + + + singlePhaseTransportModel laminarTransport(U, phi); + + autoPtr turbulence + ( + incompressible::RASModel::New(U, phi, laminarTransport) + ); + + + MRFZones mrfZones(mesh); + mrfZones.correctBoundaryVelocity(U); diff --git a/applications/solvers/incompressible/simpleSRFFoam/Make/files b/applications/solvers/incompressible/simpleSRFFoam/Make/files new file mode 100644 index 000000000..20806dff3 --- /dev/null +++ b/applications/solvers/incompressible/simpleSRFFoam/Make/files @@ -0,0 +1,3 @@ +simpleSRFFoam.C + +EXE = $(FOAM_APPBIN)/simpleSRFFoam diff --git a/applications/solvers/incompressible/simpleSRFFoam/Make/options b/applications/solvers/incompressible/simpleSRFFoam/Make/options new file mode 100644 index 000000000..02297d94c --- /dev/null +++ b/applications/solvers/incompressible/simpleSRFFoam/Make/options @@ -0,0 +1,11 @@ +EXE_INC = \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/turbulenceModels \ + -I$(LIB_SRC)/transportModels + +EXE_LIBS = \ + -lincompressibleRASModels \ + -lincompressibleTransportModels \ + -lfiniteVolume \ + -lmeshTools \ + -llduSolvers diff --git a/applications/solvers/incompressible/simpleSRFFoam/createFields.H b/applications/solvers/incompressible/simpleSRFFoam/createFields.H new file mode 100644 index 000000000..864b1dd93 --- /dev/null +++ b/applications/solvers/incompressible/simpleSRFFoam/createFields.H @@ -0,0 +1,74 @@ + Info << "Reading field p\n" << endl; + volScalarField p + ( + IOobject + ( + "p", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + Info<< "Reading field Urel\n" << endl; + volVectorField Urel + ( + IOobject + ( + "Urel", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + Info<< "Reading/calculating face flux field phi\n" << endl; + surfaceScalarField phi + ( + IOobject + ( + "phi", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + linearInterpolate(Urel) & mesh.Sf() + ); + + label pRefCell = 0; + scalar pRefValue = 0.0; + setRefCell(p, mesh.solutionDict().subDict("SIMPLE"), pRefCell, pRefValue); + + singlePhaseTransportModel laminarTransport(Urel, phi); + + autoPtr turbulence + ( + incompressible::RASModel::New(Urel, phi, laminarTransport) + ); + + Info<< "Creating SRF model\n" << endl; + autoPtr SRF + ( + SRF::SRFModel::New(Urel) + ); + + // Create Uabs as a permanent field to make it available for on-the-fly + // post-processing operations + volVectorField Uabs + ( + IOobject + ( + "Uabs", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + Urel + SRF->U() + ); + diff --git a/applications/solvers/incompressible/simpleSRFFoam/simpleSRFFoam.C b/applications/solvers/incompressible/simpleSRFFoam/simpleSRFFoam.C new file mode 100644 index 000000000..4e8768b30 --- /dev/null +++ b/applications/solvers/incompressible/simpleSRFFoam/simpleSRFFoam.C @@ -0,0 +1,131 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright held by original author + \\/ 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 2 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, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Application + simpleSRFFoam + +Description + Steady-state solver for incompressible, turbulent flow of non-Newtonian + fluids with single rotating frame. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "incompressible/singlePhaseTransportModel/singlePhaseTransportModel.H" +#include "incompressible/RAS/RASModel/RASModel.H" +#include "SRFModel.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + +# include "setRootCase.H" + +# include "createTime.H" +# include "createMesh.H" +# include "createFields.H" +# include "initContinuityErrs.H" + + //mesh.clearPrimitives(); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + for (runTime++; !runTime.end(); runTime++) + { + Info<< "Time = " << runTime.timeName() << nl << endl; + +# include "readSIMPLEControls.H" + + p.storePrevIter(); + + // Pressure-velocity SIMPLE corrector + { + // Momentum predictor + tmp UrelEqn + ( + fvm::div(phi, Urel) + + turbulence->divDevReff(Urel) + + SRF->Su() + ); + + UrelEqn().relax(); + + solve(UrelEqn() == -fvc::grad(p)); + + p.boundaryField().updateCoeffs(); + volScalarField AUrel = UrelEqn().A(); + Urel = UrelEqn().H()/AUrel; + UrelEqn.clear(); + phi = fvc::interpolate(Urel) & mesh.Sf(); + adjustPhi(phi, Urel, p); + + // Non-orthogonal pressure corrector loop + for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + { + fvScalarMatrix pEqn + ( + fvm::laplacian(1.0/AUrel, p) == fvc::div(phi) + ); + + pEqn.setReference(pRefCell, pRefValue); + pEqn.solve(); + + if (nonOrth == nNonOrthCorr) + { + phi -= pEqn.flux(); + } + } + +# include "continuityErrs.H" + + // Explicitly relax pressure for momentum corrector + p.relax(); + + // Momentum corrector + Urel -= fvc::grad(p)/AUrel; + Urel.correctBoundaryConditions(); + } + + turbulence->correct(); + + // Recalculate Uabs + Uabs = Urel + SRF->U(); + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + } + + Info<< "End\n" << endl; + + return(0); +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/MRFInterFoam/MRFInterFoam.C b/applications/solvers/multiphase/MRFInterFoam/MRFInterFoam.C new file mode 100644 index 000000000..e9415bf3e --- /dev/null +++ b/applications/solvers/multiphase/MRFInterFoam/MRFInterFoam.C @@ -0,0 +1,109 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright held by original author + \\/ 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 2 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, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Application + MRFInterFoam + +Description + Solver for 2 incompressible, isothermal immiscible fluids using a VOF + (volume of fluid) phase-fraction based interface capturing approach. + The momentum and other fluid properties are of the "mixture" and a single + momentum equation is solved. + + Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected. + + For a two-fluid approach see twoPhaseEulerFoam. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "MULES.H" +#include "subCycle.H" +#include "interfaceProperties.H" +#include "twoPhaseMixture.H" +#include "turbulenceModel.H" +#include "MRFZones.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + #include "setRootCase.H" + #include "createTime.H" + #include "createMesh.H" + #include "readGravitationalAcceleration.H" + #include "readPISOControls.H" + #include "initContinuityErrs.H" + #include "createFields.H" + #include "createMRFZones.H" + #include "readTimeControls.H" + #include "correctPhi.H" + #include "CourantNo.H" + #include "setInitialDeltaT.H" + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + while (runTime.run()) + { + #include "readPISOControls.H" + #include "readTimeControls.H" + #include "CourantNo.H" + #include "setDeltaT.H" + + runTime++; + + Info<< "Time = " << runTime.timeName() << nl << endl; + + twoPhaseProperties.correct(); + + #include "alphaEqnSubCycle.H" + + #include "UEqn.H" + + // --- PISO loop + for (int corr=0; corrcorrect(); + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + } + + Info<< "End\n" << endl; + + return 0; +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/MRFInterFoam/Make/files b/applications/solvers/multiphase/MRFInterFoam/Make/files new file mode 100644 index 000000000..9610e63ee --- /dev/null +++ b/applications/solvers/multiphase/MRFInterFoam/Make/files @@ -0,0 +1,3 @@ +MRFInterFoam.C + +EXE = $(FOAM_APPBIN)/MRFInterFoam diff --git a/applications/solvers/multiphase/MRFInterFoam/Make/options b/applications/solvers/multiphase/MRFInterFoam/Make/options new file mode 100644 index 000000000..1d614c708 --- /dev/null +++ b/applications/solvers/multiphase/MRFInterFoam/Make/options @@ -0,0 +1,15 @@ +EXE_INC = \ + -I$(FOAM_SOLVERS)/multiphase/interFoam \ + -I$(LIB_SRC)/transportModels \ + -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ + -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \ + -I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \ + -I$(LIB_SRC)/finiteVolume/lnInclude + +EXE_LIBS = \ + -linterfaceProperties \ + -lincompressibleTransportModels \ + -lincompressibleTurbulenceModel \ + -lincompressibleRASModels \ + -lincompressibleLESModels \ + -lfiniteVolume diff --git a/applications/solvers/multiphase/MRFInterFoam/UEqn.H b/applications/solvers/multiphase/MRFInterFoam/UEqn.H new file mode 100644 index 000000000..0a65664f7 --- /dev/null +++ b/applications/solvers/multiphase/MRFInterFoam/UEqn.H @@ -0,0 +1,35 @@ + surfaceScalarField muEff + ( + "muEff", + twoPhaseProperties.muf() + + fvc::interpolate(rho*turbulence->nut()) + ); + + fvVectorMatrix UEqn + ( + fvm::ddt(rho, U) + + fvm::div(rhoPhi, U) + - fvm::laplacian(muEff, U) + - (fvc::grad(U) & fvc::grad(muEff)) + //- fvc::div(muEff*(fvc::interpolate(dev(fvc::grad(U))) & mesh.Sf())) + ); + mrfZones.addCoriolis(rho, UEqn); + + UEqn.relax(); + + if (momentumPredictor) + { + solve + ( + UEqn + == + fvc::reconstruct + ( + fvc::interpolate(rho)*(g & mesh.Sf()) + + ( + fvc::interpolate(interface.sigmaK())*fvc::snGrad(alpha1) + - fvc::snGrad(pd) + ) * mesh.magSf() + ) + ); + } diff --git a/applications/solvers/multiphase/MRFInterFoam/createMRFZones.H b/applications/solvers/multiphase/MRFInterFoam/createMRFZones.H new file mode 100644 index 000000000..161446a8e --- /dev/null +++ b/applications/solvers/multiphase/MRFInterFoam/createMRFZones.H @@ -0,0 +1,2 @@ + MRFZones mrfZones(mesh); + mrfZones.correctBoundaryVelocity(U); diff --git a/applications/solvers/multiphase/MRFInterFoam/pEqn.H b/applications/solvers/multiphase/MRFInterFoam/pEqn.H new file mode 100644 index 000000000..b52c7441c --- /dev/null +++ b/applications/solvers/multiphase/MRFInterFoam/pEqn.H @@ -0,0 +1,49 @@ +{ + volScalarField rUA = 1.0/UEqn.A(); + surfaceScalarField rUAf = fvc::interpolate(rUA); + + U = rUA*UEqn.H(); + + surfaceScalarField phiU + ( + "phiU", + (fvc::interpolate(U) & mesh.Sf()) + //+ fvc::ddtPhiCorr(rUA, rho, U, phi) + ); + mrfZones.relativeFlux(phiU); + + phi = phiU + + ( + fvc::interpolate(interface.sigmaK())* + fvc::snGrad(alpha1)*mesh.magSf() + + fvc::interpolate(rho)*(g & mesh.Sf()) + )*rUAf; + adjustPhi(phi, U, p); + + for(int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + { + fvScalarMatrix pdEqn + ( + fvm::laplacian(rUAf, pd) == fvc::div(phi) + ); + + pdEqn.setReference(pdRefCell, pdRefValue); + + if (corr == nCorr-1 && nonOrth == nNonOrthCorr) + { + pdEqn.solve(mesh.solutionDict().solver(pd.name() + "Final")); + } + else + { + pdEqn.solve(mesh.solutionDict().solver(pd.name())); + } + + if (nonOrth == nNonOrthCorr) + { + phi -= pdEqn.flux(); + } + } + + U += rUA*fvc::reconstruct((phi - phiU)/rUAf); + U.correctBoundaryConditions(); +} -- 2.11.4.GIT