Post-processing solvers are uploaded.

This commit is contained in:
dhkim 2018-01-23 23:59:03 +09:00
parent ad1db37d45
commit fa2fa99978
17 changed files with 1602 additions and 0 deletions

View file

@ -0,0 +1,4 @@
NOxFoam.C
EXE = $(FOAM_USER_APPBIN)/NOxFoam_thermalNOx

View file

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I${LIB_SRC}/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
EXE_LIBS = \
-lfiniteVolume \
-lfvOptions \
-lmeshTools \
-lsampling

View file

@ -0,0 +1,97 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2016 OpenFOAM Foundation
\\/ 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 3 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, see <http://www.gnu.org/licenses/>.
Application
NOxFoam_thermalNOx
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "fvOptions.H"
#include "simpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createFields.H"
#include "createFvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
if(instantaneousRadicals==false)
{
Info<< "Partial Equilibrium Approach is selected for O, OH radicals" << endl;
}
else
{
Info<< "Instantaneous massfraction field will be used for O, OH radicals" << endl;
}
while (simple.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
#include "SnoCalc.H"
mut=rho*nut;
while (simple.correctNonOrthogonal())
{
fvScalarMatrix NOEqn
(
fvm::ddt(rho, NO)
+ fvm::div(phi, NO)
- fvm::laplacian(mut, NO)
==
Sno
);
NOEqn.relax();
fvOptions.constrain(NOEqn);
NOEqn.solve();
fvOptions.correct(NO);
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,44 @@
const scalar Wno = 0.030006; //Molecular weight of NO [kg/mol]
const scalar Wn2 = 0.028013; // N2 [kg/mol]
const scalar Wo2 = 0.031998; // O2 [kg/mol]
const scalar Wh2o = 0.018015; // H2O[kg/mol]
const scalar Wo = 0.015999; // O [kg/mol]
const scalar Woh = 0.017008; // OH [kg/mol]
forAll (mesh.cells(),celli)
{
Xno[celli]=rho[celli]*NO[celli]/Wno; //Molar concentration of NO [mol/m^3]
Xn2[celli]=rho[celli]*N2[celli]/Wn2; // N2 [mol/m^3]
Xo2[celli]=rho[celli]*O2[celli]/Wo2; // O2 [mol/m^3]
Xh2o[celli]=rho[celli]*H2O[celli]/Wh2o; // H2O[mol/m^3]
// O,OH molar concentrations
if(instantaneousRadicals==false)
{
Xo[celli]=36.64*Foam::pow(T[celli],0.5)*Foam::pow(Xo2[celli],0.5)*Foam::exp(-27123/T[celli]);
//Molar concentration of O2 [mol/m^3], Partial Equilibrium Approach
Xoh[celli]=2.129E+02*Foam::pow(T[celli],-0.57)*Foam::exp(-4595/T[celli])*Foam::pow(Xo[celli],0.5)*Foam::pow(Xh2o[celli],0.5);
//Molar concentration of OH [mol/m^3], Partial Equilibrium Approach
}
else
{
Xo[celli]=rho[celli]*O[celli]/Wo; //Molar concentration of O [mol/m^3]
Xoh[celli]=rho[celli]*OH[celli]/Woh; // OH[mol/m^3]
}
//Reaction constants
kf1[celli]=1.8E+08*Foam::exp(-38370/T[celli]);
kf2[celli]=1.8E+04*T[celli]*Foam::exp(-4680/T[celli]);
kf3[celli]=7.1E+07*Foam::exp(-450/T[celli]);
kr1[celli]=3.8E+07*Foam::exp(-425/T[celli]);
kr2[celli]=3.81E+03*T[celli]*Foam::exp(-20820/T[celli]);
Sno[celli] = Wno*2*kf1[celli]*Xo[celli]*Xn2[celli]*
((1-(kr1[celli]*kr2[celli]*Xno[celli]*Xno[celli])/(kf1[celli]*Xn2[celli]*kf2[celli]*Xo2[celli]))/
(1+(kr1[celli]*Xno[celli])/(kf2[celli]*Xo2[celli]+kf3[celli]*Xoh[celli])))/runTime.time().deltaTValue();
}

View file

@ -0,0 +1,340 @@
volScalarField T
(
IOobject
(
"T",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
Info<< "\nReading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
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(rho*U) & mesh.Sf()
);
volScalarField nut
(
IOobject
(
"nut",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField mut
(
IOobject
(
"mut",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("mut",dimensionSet(1,-1,-1,0,0,0,0),0.0)
);
volScalarField Xno
(
IOobject
(
"Xno",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xno",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xn2
(
IOobject
(
"Xn2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xn2",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xo2
(
IOobject
(
"Xo2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xo2",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xh2o
(
IOobject
(
"Xh2o",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xh2o",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xo
(
IOobject
(
"Xo",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xo",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xoh
(
IOobject
(
"Xoh",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xoh",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField NO
(
IOobject
(
"NO",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField N2
(
IOobject
(
"N2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField O2
(
IOobject
(
"O2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField H2O
(
IOobject
(
"H2O",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField O
(
IOobject
(
"O",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField OH
(
IOobject
(
"OH",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField Sno
(
IOobject
(
"Sno",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Sno",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField kf1
(
IOobject
(
"kf1",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kf1",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kf2
(
IOobject
(
"kf2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kf2",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kf3
(
IOobject
(
"kf3",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kf3",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kr1
(
IOobject
(
"kr1",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kr1",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kr2
(
IOobject
(
"kr2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kr2",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
IOdictionary modelParameter
(
IOobject
(
"modelParameter",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);
bool instantaneousRadicals(modelParameter.lookupOrDefault("instantaneousRadicals", true));

View file

@ -0,0 +1,5 @@
NOxFoam.C
EXE = $(FOAM_APPBIN)/NOxFoam_thermFuel_2

View file

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I${LIB_SRC}/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-lfvOptions \
-lsampling

View file

@ -0,0 +1,17 @@
// NH3 transport equation
while (simple.correctNonOrthogonal())
{
fvScalarMatrix NH3Eqn
(
fvm::ddt(rho, NH3)
+ fvm::div(phi, NH3)
- fvm::laplacian(mut, NH3)
==
Snh3
);
NH3Eqn.relax();
fvOptions.constrain(NH3Eqn);
NH3Eqn.solve();
fvOptions.correct(NH3);
}

View file

@ -0,0 +1,17 @@
// NO transport equation
while (simple.correctNonOrthogonal())
{
fvScalarMatrix NOEqn
(
fvm::ddt(rho, NO)
+ fvm::div(phi, NO)
- fvm::laplacian(mut, NO)
==
Sno
);
NOEqn.relax();
fvOptions.constrain(NOEqn);
NOEqn.solve();
fvOptions.correct(NO);
}

View file

@ -0,0 +1,82 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2016 OpenFOAM Foundation
\\/ 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 3 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, see <http://www.gnu.org/licenses/>.
Application
simpleReactingParcelFoam
Description
Steady state solver for compressible, turbulent flow with reacting,
multiphase particle clouds and optional sources/constraints.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "fvOptions.H"
#include "simpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createFields.H"
#include "createFvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
if(instantaneousRadicals==false)
{
Info<< "Partial Equilibrium Approach is selected for O, OH radicals" << endl;
}
else if(instantaneousRadicals==true)
{
Info<< "Instantaneous massfraction field will be used for O, OH radicals" << endl;
}
while (simple.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
mut=rho*nut;
#include "SourceCalc.H"
#include "NH3Eqn.H"
#include "NOEqn.H"
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,94 @@
const scalar Wno = 0.030006; //Molecular weight of NO [kg/mol]
const scalar Wn2 = 0.028013; // N2 [kg/mol]
const scalar Wo2 = 0.031998; // O2 [kg/mol]
const scalar Wh2o = 0.018015; // H2O[kg/mol]
const scalar Wo = 0.015999; // O [kg/mol]
const scalar Woh = 0.017008; // OH [kg/mol]
const scalar Wnh3 = 0.017031; // NH3[kg/mol]
const scalar Wn = 0.014007; // N [kg/mol]
const scalar RR = 8.314; //Universal gas constant [J/(mol*K)]
const scalar A1 = 4.0E+06; // [1/s]
const scalar A2 = 1.8E+08; // [1/s]
const scalar E1 = 133947.2; // [J/mol]
const scalar E2 = 113017.95; // [J/mol]
forAll (mesh.cells(),celli)
{
// molar concentrations
Xno[celli]=rho[celli]*NO[celli]/Wno; //Molar concentration of NO [mol/m^3]
Xn2[celli]=rho[celli]*N2[celli]/Wn2; // N2 [mol/m^3]
Xo2[celli]=rho[celli]*O2[celli]/Wo2; // O2 [mol/m^3]
Xh2o[celli]=rho[celli]*H2O[celli]/Wh2o; // H2O[mol/m^3]
Xnh3[celli]=rho[celli]*NH3[celli]/Wnh3; // nh3[mol/m^3]
// O,OH molar concentrations
if(instantaneousRadicals==false)
{
Xo[celli]=36.64*Foam::pow(T[celli],0.5)*Foam::pow(Xo2[celli],0.5)*Foam::exp(-27123/T[celli]);
//Molar concentration of O2 [mol/m^3], Partial Equilibrium Approach
Xoh[celli]=2.129E+02*Foam::pow(T[celli],-0.57)*Foam::exp(-4595/T[celli])*Foam::pow(Xo[celli],0.5)*Foam::pow(Xh2o[celli],0.5);
//Molar concentration of OH [mol/m^3], Partial Equilibrium Approach
}
else
{
Xo[celli]=rho[celli]*O[celli]/Wo; //Molar concentration of O [mol/m^3]
Xoh[celli]=rho[celli]*OH[celli]/Woh; // OH[mol/m^3]
}
// NH3
//Oxygen Reaction Order, a
if (Xo2[celli] <= 4.1E-03) {
a[celli] = 1.0;
}
else if ( Xo2[celli] > 4.1E-03 and Xo2[celli] <= 1.11E-02 ) {
a[celli] = -3.95-0.9*Foam::log(Xo2[celli]);
}
else if ( Xo2[celli] > 1.11-02 and Xo2[celli] <= 0.03 ) {
a[celli] = -0.35-0.1*Foam::log(Xo2[celli]);
}
else if ( Xo2[celli] > 0.03 ) {
a[celli] = 0;
}
//Conversion rate of NH3
R1[celli] = A1*Xnh3[celli]*Foam::pow(Xo2[celli],a[celli])*Foam::exp(-E1/(RR*T[celli])); // [1/s]
R2[celli] = A2*Xnh3[celli]*Xno[celli]*Foam::exp(-E2/(RR*T[celli])); // [1/s]
//Source and sink of NH3
Snh3_p[celli] = (Sfuel1[celli]+Sfuel2[celli])*Yn[celli]*Wnh3/Wn/mesh.V()[celli]; //NH3 production
Snh3_1[celli] = -R1[celli]*Wnh3*p[celli]/(RR*T[celli]); //NH3 consumption -> NO (oxidation)
Snh3_2[celli] = -R2[celli]*Wnh3*p[celli]/(RR*T[celli]); //NH3 consumption -> N2 (reduction)
//Sum of Sources (NH3 production, consumption 1, and consumption 2 )
Snh3[celli] = Snh3_p[celli] +Snh3_1[celli] +Snh3_2[celli];
//Thermal NO
//Reaction rate constants [m^3/(mol*s)]
kf1[celli]=1.8E+08*Foam::exp(-38370/T[celli]);
kf2[celli]=1.8E+04*T[celli]*Foam::exp(-4680/T[celli]);
kf3[celli]=7.1E+07*Foam::exp(-450/T[celli]);
kr1[celli]=3.8E+07*Foam::exp(-425/T[celli]);
kr2[celli]=3.81E+03*T[celli]*Foam::exp(-20820/T[celli]);
//Calculation of source term. (Thermal NO)
SthermNO[celli] = Wno*2*kf1[celli]*Xo[celli]*Xn2[celli]*
((1-(kr1[celli]*kr2[celli]*Xno[celli]*Xno[celli])/
(kf1[celli]*Xn2[celli]*kf2[celli]*Xo2[celli]))/
(1+(kr1[celli]*Xno[celli])/
(kf2[celli]*Xo2[celli]+kf3[celli]*Xoh[celli])))
/runTime.time().deltaTValue();
// Fuel NO
//Calculation of source term. (Fuel NO)
SfuelNO_1[celli] = R1[celli]*Wno*p[celli]/(RR*T[celli]); // NH3 + O2 -> NO (source)
SfuelNO_2[celli] =-R2[celli]*Wno*p[celli]/(RR*T[celli]); // NH3 + NO -> N2 (sink)
SfuelNO[celli] = SfuelNO_1[celli] + SfuelNO_2[celli]; // Sum of fuel NO source terms [kg/(m^3*s)]
// Sum of NO Sources (source of thermal NO and fuel NO)
Sno[celli] = SthermNO[celli] + SfuelNO[celli]; // [kg/(m^3*s)]
}

View file

@ -0,0 +1,591 @@
volScalarField T
(
IOobject
(
"T",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField p
(
IOobject
(
"p",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
surfaceScalarField phi
(
IOobject
(
"phi",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
linearInterpolate(rho*U) & mesh.Sf()
);
volScalarField nut
(
IOobject
(
"nut",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField mut
(
IOobject
(
"mut",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("mut",dimensionSet(1,-1,-1,0,0,0,0),0.0)
);
volScalarField Xno
(
IOobject
(
"Xno",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xno",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xn2
(
IOobject
(
"Xn2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xn2",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xo2
(
IOobject
(
"Xo2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xo2",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xh2o
(
IOobject
(
"Xh2o",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xh2o",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xnh3
(
IOobject
(
"Xnh3",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xnh3",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xo
(
IOobject
(
"Xo",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xo",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField Xoh
(
IOobject
(
"Xoh",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("Xoh",dimensionSet(0,-3,0,0,1,0,0),0.0)
);
volScalarField NO
(
IOobject
(
"NO",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField N2
(
IOobject
(
"N2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField O2
(
IOobject
(
"O2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField H2O
(
IOobject
(
"H2O",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField O
(
IOobject
(
"O",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField OH
(
IOobject
(
"OH",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField NH3
(
IOobject
(
"NH3",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField Sno
(
IOobject
(
"Sno",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Sno",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField SthermNO
(
IOobject
(
"SthermNO",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("SthermNO",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField SfuelNO
(
IOobject
(
"SfuelNO",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("SfuelNO",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField SfuelNO_1
(
IOobject
(
"SfuelNO_1",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("SfuelNO_1",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField SfuelNO_2
(
IOobject
(
"SfuelNO_2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("SfuelNO_2",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField a
(
IOobject
(
"a",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("a",dimensionSet(0,0,0,0,0,0,0),0.0)
);
volScalarField R1
(
IOobject
(
"R1",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("R1",dimensionSet(0,0,-1,0,0,0,0),0.0)
);
volScalarField R2
(
IOobject
(
"R2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("R1",dimensionSet(0,0,-1,0,0,0,0),0.0)
);
volScalarField Snh3_p
(
IOobject
(
"Snh3_p",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Snh3_p",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField Snh3_1
(
IOobject
(
"Snh3_1",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Snh3_1",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField Snh3_2
(
IOobject
(
"Snh3_2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Snh3_2",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField Snh3
(
IOobject
(
"Snh3",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Snh3",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField kf1
(
IOobject
(
"kf1",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kf1",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kf2
(
IOobject
(
"kf2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kf2",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kf3
(
IOobject
(
"kf3",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kf3",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kr1
(
IOobject
(
"kr1",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kr1",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField kr2
(
IOobject
(
"kr2",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("kr2",dimensionSet(0,3,-1,0,-1,0,0),0.0)
);
volScalarField Sfuel1
(
IOobject
(
"Sfuel1",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField Sfuel2
(
IOobject
(
"Sfuel2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
IOdictionary fuelNitrogenMassFraction
(
IOobject
(
"fuelNitrogenMassFraction",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);
volScalarField Yn
(
IOobject
(
"Yn",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensionedScalar(fuelNitrogenMassFraction.lookup("Yn"))
);
IOdictionary modelParameter
(
IOobject
(
"modelParameter",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);
bool instantaneousRadicals(modelParameter.lookupOrDefault("instantaneousRadicals", true));

View file

@ -0,0 +1,4 @@
SOxFoam.C
EXE = $(FOAM_USER_APPBIN)/SOxFoam

View file

@ -0,0 +1,11 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I${LIB_SRC}/meshTools/lnInclude \
-I$(LIB_SRC)/sampling/lnInclude \
EXE_LIBS = \
-lfiniteVolume \
-lmeshTools \
-lfvOptions \
-lsampling

View file

@ -0,0 +1,89 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2016 OpenFOAM Foundation
\\/ 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 3 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, see <http://www.gnu.org/licenses/>.
Application
simpleReactingParcelFoam
Description
Steady state solver for compressible, turbulent flow with reacting,
multiphase particle clouds and optional sources/constraints.
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "fvOptions.H"
#include "simpleControl.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
#include "createControl.H"
#include "createFields.H"
#include "createFvOptions.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (simple.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
#include "Sso2Calc.H"
mut=rho*nut;
while (simple.correctNonOrthogonal())
{
fvScalarMatrix SO2Eqn
(
fvm::ddt(rho, SO2)
+ fvm::div(phi, SO2)
- fvm::laplacian(mut, SO2)
==
Sso2
);
SO2Eqn.relax();
fvOptions.constrain(SO2Eqn);
SO2Eqn.solve();
fvOptions.correct(SO2);
}
runTime.write();
Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View file

@ -0,0 +1,11 @@
const scalar Wso2 = 0.064066; //Molecular weight of SO2[kg/mol]
const scalar Ws = 0.032065; // S [kg/mol]
forAll (mesh.cells(),celli)
{
Sso2[celli] = (Sfuel1[celli]+Sfuel2[celli])*Ys[celli]*Wso2/Ws/mesh.V()[celli];
}

View file

@ -0,0 +1,176 @@
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
//thermo.rho()
mesh
);
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
surfaceScalarField phi
(
IOobject
(
"phi",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
linearInterpolate(rho*U) & mesh.Sf()
);
volScalarField nut
(
IOobject
(
"nut",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField mut
(
IOobject
(
"mut",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("mut",dimensionSet(1,-1,-1,0,0,0,0),0.0)
);
//#include "compressibleCreatePhi.H"
//
//mesh.setFluxRequired(p.name());
//
//dimensionedScalar rhoMax
//(
// dimensionedScalar::lookupOrDefault
// (
// "rhoMax",
// simple.dict(),
// dimDensity,
// GREAT
// )
//);
//
//dimensionedScalar rhoMin
//(
// dimensionedScalar::lookupOrDefault
// (
// "rhoMin",
// simple.dict(),
// dimDensity,
// 0
// )
//);
volScalarField SO2
(
IOobject
(
"SO2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField Sso2
(
IOobject
(
"Sso2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Sso2",dimensionSet(1,-3,-1,0,0,0,0),0.0)
);
volScalarField Sfuel1
(
IOobject
(
"Sfuel1",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField Sfuel2
(
IOobject
(
"Sfuel2",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
IOdictionary fuelSulfurMassFraction
(
IOobject
(
"fuelSulfurMassFraction",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);
volScalarField Ys
(
IOobject
(
"Ys",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar(fuelSulfurMassFraction.lookup("Ys"))
);