From fd71f367d7dfbb29dec71a9a557dc12b313f7b09 Mon Sep 17 00:00:00 2001 From: asanaz Date: Thu, 6 Sep 2018 12:54:35 +0200 Subject: [PATCH] rCFD solver for turbulent single-phase transport --- .../solvers/recSolverTurbTransport/Make/files | 3 + .../recSolverTurbTransport/Make/options | 27 +++ .../solvers/recSolverTurbTransport/TEq.H | 17 ++ .../recSolverTurbTransport/createFields.H | 174 ++++++++++++++++++ .../recSolverTurbTransport/readFields.H | 14 ++ .../recSolverTurbTransport.C | 113 ++++++++++++ 6 files changed, 348 insertions(+) create mode 100755 applications/solvers/recSolverTurbTransport/Make/files create mode 100755 applications/solvers/recSolverTurbTransport/Make/options create mode 100755 applications/solvers/recSolverTurbTransport/TEq.H create mode 100755 applications/solvers/recSolverTurbTransport/createFields.H create mode 100755 applications/solvers/recSolverTurbTransport/readFields.H create mode 100755 applications/solvers/recSolverTurbTransport/recSolverTurbTransport.C diff --git a/applications/solvers/recSolverTurbTransport/Make/files b/applications/solvers/recSolverTurbTransport/Make/files new file mode 100755 index 00000000..c984de5c --- /dev/null +++ b/applications/solvers/recSolverTurbTransport/Make/files @@ -0,0 +1,3 @@ +recSolverTurbTransport.C + +EXE=$(CFDEM_APP_DIR)/recSolverTurbTransport diff --git a/applications/solvers/recSolverTurbTransport/Make/options b/applications/solvers/recSolverTurbTransport/Make/options new file mode 100755 index 00000000..686f676c --- /dev/null +++ b/applications/solvers/recSolverTurbTransport/Make/options @@ -0,0 +1,27 @@ +include $(CFDEM_ADD_LIBS_DIR)/additionalLibs + +EXE_INC = \ + -I$(CFDEM_OFVERSION_DIR) \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \ + -I$(LIB_SRC)/transportModels \ + -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \ + -I$(CFDEM_SRC_DIR)/lagrangian/cfdemParticle/lnInclude \ + -I$(CFDEM_SRC_DIR)/lagrangian/cfdemParticle/cfdTools \ + -I$(CFDEM_SRC_DIR)/recurrence/lnInclude \ + -I$(CFDEM_SRC_DIR)/lagrangian/cfdemParticle/derived/cfdemCloudRec \ + +EXE_LIBS = \ + -L$(CFDEM_LIB_DIR)\ + -lrecurrence \ + -lturbulenceModels \ + -lincompressibleTurbulenceModels \ + -lincompressibleTransportModels \ + -lfiniteVolume \ + -lmeshTools \ + -lfvOptions \ + -l$(CFDEM_LIB_NAME) \ + $(CFDEM_ADD_LIB_PATHS) \ + $(CFDEM_ADD_LIBS) diff --git a/applications/solvers/recSolverTurbTransport/TEq.H b/applications/solvers/recSolverTurbTransport/TEq.H new file mode 100755 index 00000000..e28a172d --- /dev/null +++ b/applications/solvers/recSolverTurbTransport/TEq.H @@ -0,0 +1,17 @@ + +volScalarField alphaEff("alphaEff", turbulence->nu()/Sc + dU2/Sct); + +TEqn = +( + fvm::ddt(T) + + fvm::div(phiRec, T) + - fvm::laplacian(alphaEff, T) + == + fvOptions(T) +); + +TEqn.relax(relaxCoeff); + +fvOptions.constrain(TEqn); + +TEqn.solve(); diff --git a/applications/solvers/recSolverTurbTransport/createFields.H b/applications/solvers/recSolverTurbTransport/createFields.H new file mode 100755 index 00000000..3b9b3af1 --- /dev/null +++ b/applications/solvers/recSolverTurbTransport/createFields.H @@ -0,0 +1,174 @@ + // dummy fields + Info<< "\nCreating dummy pressure and density fields\n" << endl; + volScalarField p + ( + IOobject + ( + "p", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh, + dimensionedScalar("p", dimensionSet(1, 2, -2, 0, 0), 1.0) + ); + + volScalarField rho + ( + IOobject + ( + "rho", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh, + dimensionedScalar("rho", dimensionSet(1, -3, 0, 0, 0), 1.0) + ); + + // recurrence fields + Info<< "\nCreating recurrence fields.\n" << endl; + + volVectorField URec + ( + IOobject + ( + "URec", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + volScalarField U2Rec + ( + IOobject + ( + "U2Rec", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + // calculated fields + Info<< "\nCreating fields subject to calculation\n" << endl; + + volScalarField delta + ( + IOobject + ( + "delta", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + mesh, + dimensionedScalar("delta", dimLength, 0.0) + ); + + delta.primitiveFieldRef()=pow(mesh.V(),1.0/3.0); + delta.write(); + + + Info<< "\ncreating dU2\n" << endl; + + volScalarField dU2 + ( + IOobject + ( + "dU2", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + sqrt(0.5*mag(U2Rec - magSqr(URec)))*delta*0.094 + ); + + forAll(dU2, cellI) + { + if (U2Rec[cellI]-magSqr(URec[cellI]) < 0.0) + { + dU2[cellI] = 0.0; + } + } + + dU2.write(); + +Info<< "Calculating face flux field phiRec\n" << endl; +surfaceScalarField phiRec + ( + IOobject + ( + "phiRec", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + linearInterpolate(URec) & mesh.Sf() + ); + + phiRec.write(); + + singlePhaseTransportModel laminarTransport(URec, phiRec); + + autoPtr turbulence + ( + incompressible::turbulenceModel::New(URec, phiRec, laminarTransport) + ); + + dimensionedScalar Sc("Sc", dimless, laminarTransport); + dimensionedScalar Sct("Sct", dimless, laminarTransport); + + // create concentration field + Info<< "Creating scalar transport field\n" << endl; + + volScalarField T + ( + IOobject + ( + "T", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + fvScalarMatrix TEqn(T, dimless*dimVolume/(dimTime)); + + scalar relaxCoeff(0.0); + + Info<< "reading clockProperties\n" << endl; + + IOdictionary clockProperties + ( + IOobject + ( + "clockProperties", + mesh.time().constant(), + mesh, + IOobject::MUST_READ, + IOobject::NO_WRITE + ) + ); + + autoPtr myClock + ( + clockModel::New + ( + clockProperties, + mesh.time() + ) + ); + diff --git a/applications/solvers/recSolverTurbTransport/readFields.H b/applications/solvers/recSolverTurbTransport/readFields.H new file mode 100755 index 00000000..13dad011 --- /dev/null +++ b/applications/solvers/recSolverTurbTransport/readFields.H @@ -0,0 +1,14 @@ + +recurrenceBase.recM().exportVolScalarField("U2Mean",U2Rec); +recurrenceBase.recM().exportVolVectorField("UMean",URec); +phiRec=linearInterpolate(URec) & mesh.Sf(); + +dU2=sqrt(0.5*mag(U2Rec - magSqr(URec)))*delta*0.094; + +forAll(dU2, cellI) +{ + if (U2Rec[cellI]-magSqr(URec[cellI]) < 0.0) + { + dU2[cellI] = 0.0; + } +} diff --git a/applications/solvers/recSolverTurbTransport/recSolverTurbTransport.C b/applications/solvers/recSolverTurbTransport/recSolverTurbTransport.C new file mode 100755 index 00000000..e104b1b6 --- /dev/null +++ b/applications/solvers/recSolverTurbTransport/recSolverTurbTransport.C @@ -0,0 +1,113 @@ +/*---------------------------------------------------------------------------*\ + CFDEMcoupling academic - Open Source CFD-DEM coupling + + Contributing authors: + Thomas Lichtenegger, Gerhard Holzinger + Copyright (C) 2015- Johannes Kepler University, Linz +------------------------------------------------------------------------------- +License + This file is part of CFDEMcoupling academic. + + CFDEMcoupling academic 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 3 of the License, or + (at your option) any later version. + + CFDEMcoupling academic 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 CFDEMcoupling academic. If not, see . + +Application + Turbulent Transport Solver Recurrence + +Description + Solves a transport equation for a passive scalar on a single-phase solution + for a solver based on recurrence statistics + +Rules + Solution data to compute the recurrence statistics from, needs to + reside in $CASE_ROOT/dataBase + Time step data in dataBase needs to be evenly spaced in time + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "singlePhaseTransportModel.H" +#include "turbulentTransportModel.H" +#include "fvOptions.H" + +#include "recBase.H" +#include "recModel.H" + +#include "clockModel.H" +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + #include "postProcess.H" + #include "setRootCase.H" + #include "createTime.H" + #include "createMesh.H" + #include "createControl.H" + #include "createFields.H" + #include "createFvOptions.H" + + recBase recurrenceBase(mesh); + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nCalculating particle trajectories based on recurrence statistics\n" << endl; + + label recTimeIndex(0); + scalar recTimeStep_=recurrenceBase.recM().recTimeStep(); + + while (runTime.run()) + { + + myClock().start(1,"Global"); + + runTime++; + + Info<< "Time = " << runTime.timeName() << nl << endl; + + myClock().start(2,"fieldUpdate"); + + if ( runTime.timeOutputValue() - (recTimeIndex+1)*recTimeStep_ + 1.0e-5 > 0.0 ) + { + Info << "Updating fields at run time " << runTime.timeOutputValue() + << " corresponding to recurrence time " << (recTimeIndex+1)*recTimeStep_ << ".\n" << endl; + recurrenceBase.updateRecFields(); + #include "readFields.H" + recTimeIndex++; + } + + myClock().stop("fieldUpdate"); + + myClock().start(3,"speciesEqn"); + #include "TEq.H" + myClock().stop("speciesEqn"); + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + + myClock().stop("Global"); + + } + + myClock().evalPar(); + myClock().normHist(); + + Info<< "End\n" << endl; + + return 0; +} + + +// ************************************************************************* //