mirror of
https://develop.openfoam.com/Development/openfoam.git
synced 2025-11-28 03:28:01 +00:00
ENH: Optimization of
1) fvDOM BC. greyDiffusiveRadiationMixedFvPatchScalarField.C 2) radiationProperties and fvSolution for rays in fvDOM model 3) radiativeIntensityRay.C using weighted residuals.
This commit is contained in:
@ -44,9 +44,7 @@ bool Foam::pairPatchAgglomeration::continueAgglomerating
|
||||
{
|
||||
// Check the need for further agglomeration on all processors
|
||||
label localnCoarseFaces = nCoarseFaces;
|
||||
// reduce(localnCoarseFaces, sumOp<label>());
|
||||
bool contAgg = localnCoarseFaces >= nFacesInCoarsestLevel_;
|
||||
//reduce(contAgg, andOp<bool>());
|
||||
return contAgg;
|
||||
}
|
||||
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -189,51 +189,40 @@ updateCoeffs()
|
||||
}
|
||||
|
||||
scalarField& Iw = *this;
|
||||
const vectorField n(patch().Sf()/patch().magSf());
|
||||
const vectorField n(patch().nf());
|
||||
|
||||
radiativeIntensityRay& ray =
|
||||
const_cast<radiativeIntensityRay&>(dom.IRay(rayId));
|
||||
|
||||
ray.Qr().boundaryField()[patchI] += Iw*(n & ray.dAve());
|
||||
const scalarField nAve(n & ray.dAve());
|
||||
|
||||
scalarField temissivity = emissivity();
|
||||
ray.Qr().boundaryField()[patchI] += Iw*nAve;
|
||||
|
||||
const scalarField temissivity = emissivity();
|
||||
|
||||
scalarField& Qem = ray.Qem().boundaryField()[patchI];
|
||||
scalarField& Qin = ray.Qin().boundaryField()[patchI];
|
||||
|
||||
const vector& myRayId = dom.IRay(rayId).d();
|
||||
|
||||
const scalarField& Ir = dom.Qin();
|
||||
|
||||
forAll(Iw, faceI)
|
||||
{
|
||||
scalar Ir = 0.0;
|
||||
|
||||
for (label rayI=0; rayI < dom.nRay(); rayI++)
|
||||
{
|
||||
const vector& d = dom.IRay(rayI).d();
|
||||
|
||||
const scalarField& IFace =
|
||||
dom.IRay(rayI).ILambda(lambdaId).boundaryField()[patchI];
|
||||
|
||||
if ((-n[faceI] & d) < 0.0)
|
||||
{
|
||||
// q into the wall
|
||||
const vector& dAve = dom.IRay(rayI).dAve();
|
||||
Ir += IFace[faceI]*mag(n[faceI] & dAve);
|
||||
}
|
||||
}
|
||||
|
||||
const vector& d = dom.IRay(rayId).d();
|
||||
|
||||
if ((-n[faceI] & d) > 0.0)
|
||||
if ((-n[faceI] & myRayId) > 0.0)
|
||||
{
|
||||
// direction out of the wall
|
||||
refGrad()[faceI] = 0.0;
|
||||
valueFraction()[faceI] = 1.0;
|
||||
refValue()[faceI] =
|
||||
(
|
||||
Ir*(scalar(1.0) - temissivity[faceI])
|
||||
Ir[faceI]*(scalar(1.0) - temissivity[faceI])
|
||||
+ temissivity[faceI]*physicoChemical::sigma.value()
|
||||
* pow4(Tp[faceI])
|
||||
)/pi;
|
||||
|
||||
// Emmited heat flux from this ray direction
|
||||
ray.Qem().boundaryField()[patchI][faceI] =
|
||||
refValue()[faceI]*(n[faceI] & ray.dAve());
|
||||
// Emmited heat flux from this ray direction
|
||||
Qem[faceI] = refValue()[faceI]*nAve[faceI];
|
||||
}
|
||||
else
|
||||
{
|
||||
@ -242,9 +231,8 @@ updateCoeffs()
|
||||
refGrad()[faceI] = 0.0;
|
||||
refValue()[faceI] = 0.0; //not used
|
||||
|
||||
// Incident heat flux on this ray direction
|
||||
ray.Qin().boundaryField()[patchI][faceI] =
|
||||
Iw[faceI]*(n[faceI] & ray.dAve());
|
||||
// Incident heat flux on this ray direction
|
||||
Qin[faceI] = Iw[faceI]*nAve[faceI];
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -27,6 +27,7 @@ License
|
||||
#include "absorptionEmissionModel.H"
|
||||
#include "scatterModel.H"
|
||||
#include "constants.H"
|
||||
#include "fvm.H"
|
||||
|
||||
using namespace Foam::constant;
|
||||
using namespace Foam::constant::mathematical;
|
||||
@ -172,10 +173,39 @@ void Foam::radiation::fvDOM::initialise()
|
||||
Info<< "fvDOM : Allocated " << IRay_.size()
|
||||
<< " rays with average orientation:" << nl;
|
||||
|
||||
forAll(IRay_, i)
|
||||
if (cacheDiv_)
|
||||
{
|
||||
Info<< '\t' << IRay_[i].I().name()
|
||||
<< '\t' << IRay_[i].dAve() << nl;
|
||||
Info<< "Caching div fvMatrix..."<< endl;
|
||||
for (label lambdaI = 0; lambdaI < nLambda_; lambdaI++)
|
||||
{
|
||||
fvRayDiv_[lambdaI].setSize(nRay_);
|
||||
|
||||
forAll(IRay_, rayId)
|
||||
{
|
||||
const surfaceScalarField Ji = IRay_[rayId].dAve() & mesh_.Sf();
|
||||
const volScalarField& iRayLambdaI =
|
||||
IRay_[rayId].ILambda(lambdaI);
|
||||
|
||||
fvRayDiv_[lambdaI].set
|
||||
(
|
||||
rayId,
|
||||
new fvScalarMatrix
|
||||
(
|
||||
fvm::div(Ji, iRayLambdaI, "div(Ji,Ii_h)")
|
||||
)
|
||||
);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
forAll(IRay_, rayId)
|
||||
{
|
||||
if (omegaMax_ < IRay_[rayId].omega())
|
||||
{
|
||||
omegaMax_ = IRay_[rayId].omega();
|
||||
}
|
||||
Info<< '\t' << IRay_[rayId].I().name() << " : " << "omega : "
|
||||
<< '\t' << IRay_[rayId].omega() << nl;
|
||||
}
|
||||
|
||||
Info<< endl;
|
||||
@ -194,7 +224,7 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
|
||||
"G",
|
||||
mesh_.time().timeName(),
|
||||
mesh_,
|
||||
IOobject::READ_IF_PRESENT,
|
||||
IOobject::NO_READ,
|
||||
IOobject::AUTO_WRITE
|
||||
),
|
||||
mesh_,
|
||||
@ -247,7 +277,7 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
|
||||
mesh_.time().timeName(),
|
||||
mesh_,
|
||||
IOobject::NO_READ,
|
||||
IOobject::NO_WRITE
|
||||
IOobject::AUTO_WRITE
|
||||
),
|
||||
mesh_,
|
||||
dimensionedScalar("a", dimless/dimLength, 0.0)
|
||||
@ -260,7 +290,10 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
|
||||
blackBody_(nLambda_, T),
|
||||
IRay_(0),
|
||||
convergence_(coeffs_.lookupOrDefault<scalar>("convergence", 0.0)),
|
||||
maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50))
|
||||
maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50)),
|
||||
fvRayDiv_(nLambda_),
|
||||
cacheDiv_(coeffs_.lookupOrDefault<bool>("cacheDiv", false)),
|
||||
omegaMax_(0)
|
||||
{
|
||||
initialise();
|
||||
}
|
||||
@ -346,7 +379,10 @@ Foam::radiation::fvDOM::fvDOM
|
||||
blackBody_(nLambda_, T),
|
||||
IRay_(0),
|
||||
convergence_(coeffs_.lookupOrDefault<scalar>("convergence", 0.0)),
|
||||
maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50))
|
||||
maxIter_(coeffs_.lookupOrDefault<label>("maxIter", 50)),
|
||||
fvRayDiv_(nLambda_),
|
||||
cacheDiv_(coeffs_.lookupOrDefault<bool>("cacheDiv", false)),
|
||||
omegaMax_(0)
|
||||
{
|
||||
initialise();
|
||||
}
|
||||
@ -364,7 +400,6 @@ bool Foam::radiation::fvDOM::read()
|
||||
if (radiationModel::read())
|
||||
{
|
||||
// Only reading solution parameters - not changing ray geometry
|
||||
|
||||
coeffs_.readIfPresent("convergence", convergence_);
|
||||
coeffs_.readIfPresent("maxIter", maxIter_);
|
||||
|
||||
@ -383,19 +418,30 @@ void Foam::radiation::fvDOM::calculate()
|
||||
|
||||
updateBlackBodyEmission();
|
||||
|
||||
// Set rays convergence false
|
||||
List<bool> rayIdConv(nRay_, false);
|
||||
|
||||
scalar maxResidual = 0.0;
|
||||
label radIter = 0;
|
||||
do
|
||||
{
|
||||
Info<< "Radiation solver iter: " << radIter << endl;
|
||||
|
||||
radIter++;
|
||||
maxResidual = 0.0;
|
||||
forAll(IRay_, rayI)
|
||||
{
|
||||
maxResidual = 0.0;
|
||||
scalar maxBandResidual = IRay_[rayI].correct();
|
||||
maxResidual = max(maxBandResidual, maxResidual);
|
||||
}
|
||||
if (!rayIdConv[rayI])
|
||||
{
|
||||
scalar maxBandResidual = IRay_[rayI].correct();
|
||||
maxResidual = max(maxBandResidual, maxResidual);
|
||||
|
||||
Info<< "Radiation solver iter: " << radIter << endl;
|
||||
if (maxBandResidual < convergence_)
|
||||
{
|
||||
rayIdConv[rayI] = true;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
} while (maxResidual > convergence_ && radIter < maxIter_);
|
||||
|
||||
@ -433,7 +479,7 @@ Foam::radiation::fvDOM::Ru() const
|
||||
const DimensionedField<scalar, volMesh> E =
|
||||
absorptionEmission_->ECont()().dimensionedInternalField();
|
||||
const DimensionedField<scalar, volMesh> a =
|
||||
a_.dimensionedInternalField(); //absorptionEmission_->aCont()()
|
||||
a_.dimensionedInternalField();
|
||||
|
||||
return a*G - E;
|
||||
}
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -30,6 +30,7 @@ Description
|
||||
directions in a participating media, not including scatter.
|
||||
|
||||
Available absorption models:
|
||||
constantAbsorptionEmission
|
||||
greyMeanAbsoprtionEmission
|
||||
wideBandAbsorptionEmission
|
||||
|
||||
@ -37,9 +38,15 @@ Description
|
||||
\verbatim
|
||||
fvDOMCoeffs
|
||||
{
|
||||
nPhi 1; // azimuthal angles in PI/2 on X-Y.(from Y to X)
|
||||
nTheta 2; // polar angles in PI (from Z to X-Y plane)
|
||||
convergence 1e-4; // convergence criteria for radiation iteration
|
||||
nPhi 4; // azimuthal angles in PI/2 on X-Y.
|
||||
//(from Y to X)
|
||||
nTheta 0; // polar angles in PI (from Z to X-Y plane)
|
||||
convergence 1e-3; // convergence criteria for radiation
|
||||
//iteration
|
||||
maxIter 4; // maximum number of iterations
|
||||
cacheDiv true; // cache the div of the RTE equation.
|
||||
//NOTE: Caching div is "only" accurate if the upwind scheme is used
|
||||
//in div(Ji,Ii_h)
|
||||
}
|
||||
|
||||
solverFreq 1; // Number of flow iterations per radiation iteration
|
||||
@ -61,6 +68,7 @@ SourceFiles
|
||||
|
||||
#include "radiativeIntensityRay.H"
|
||||
#include "radiationModel.H"
|
||||
#include "fvMatrices.H"
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
@ -79,6 +87,7 @@ class fvDOM
|
||||
{
|
||||
// Private data
|
||||
|
||||
|
||||
//- Incident radiation [W/m2]
|
||||
volScalarField G_;
|
||||
|
||||
@ -121,6 +130,15 @@ class fvDOM
|
||||
//- Maximum number of iterations
|
||||
scalar maxIter_;
|
||||
|
||||
//- List of cached fvMatrices for rays
|
||||
List<PtrList<fvScalarMatrix> >fvRayDiv_;
|
||||
|
||||
//- Cache convection div matrix
|
||||
bool cacheDiv_;
|
||||
|
||||
//- Maximum omega weight
|
||||
scalar omegaMax_;
|
||||
|
||||
|
||||
// Private Member Functions
|
||||
|
||||
@ -229,6 +247,19 @@ public:
|
||||
|
||||
//- Const access to black body
|
||||
inline const blackBodyEmission& blackBody() const;
|
||||
|
||||
//- Const access to cached fvMatrix
|
||||
inline const fvScalarMatrix& fvRayDiv
|
||||
(
|
||||
const label lambdaI,
|
||||
const label rayId
|
||||
) const;
|
||||
|
||||
//- Caching div(Ji, Ilamda)
|
||||
inline bool cacheDiv() const;
|
||||
|
||||
//- Return omegaMax
|
||||
inline scalar omegaMax() const;
|
||||
};
|
||||
|
||||
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -102,6 +102,7 @@ inline const Foam::volScalarField& Foam::radiation::fvDOM::Qem() const
|
||||
return Qem_;
|
||||
}
|
||||
|
||||
|
||||
inline const Foam::radiation::blackBodyEmission&
|
||||
Foam::radiation::fvDOM::blackBody() const
|
||||
{
|
||||
@ -109,4 +110,26 @@ Foam::radiation::fvDOM::blackBody() const
|
||||
}
|
||||
|
||||
|
||||
inline const Foam::fvScalarMatrix& Foam::radiation::fvDOM::fvRayDiv
|
||||
(
|
||||
const label rayId,
|
||||
const label lambdaI
|
||||
) const
|
||||
{
|
||||
return fvRayDiv_[lambdaI][rayId];
|
||||
}
|
||||
|
||||
|
||||
inline bool Foam::radiation::fvDOM::cacheDiv() const
|
||||
{
|
||||
return cacheDiv_;
|
||||
}
|
||||
|
||||
|
||||
inline Foam::scalar Foam::radiation::fvDOM::omegaMax() const
|
||||
{
|
||||
return omegaMax_;
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -113,7 +113,8 @@ Foam::radiation::radiativeIntensityRay::radiativeIntensityRay
|
||||
phi_(phi),
|
||||
omega_(0.0),
|
||||
nLambda_(nLambda),
|
||||
ILambda_(nLambda)
|
||||
ILambda_(nLambda),
|
||||
myRayId_(rayId)
|
||||
{
|
||||
scalar sinTheta = Foam::sin(theta);
|
||||
scalar cosTheta = Foam::cos(theta);
|
||||
@ -135,7 +136,6 @@ Foam::radiation::radiativeIntensityRay::radiativeIntensityRay
|
||||
0.5*deltaPhi*Foam::sin(2.0*theta)*Foam::sin(deltaTheta)
|
||||
);
|
||||
|
||||
|
||||
autoPtr<volScalarField> IDefaultPtr;
|
||||
|
||||
forAll(ILambda_, lambdaI)
|
||||
@ -213,29 +213,51 @@ Foam::scalar Foam::radiation::radiativeIntensityRay::correct()
|
||||
{
|
||||
const volScalarField& k = dom_.aLambda(lambdaI);
|
||||
|
||||
const surfaceScalarField Ji(dAve_ & mesh_.Sf());
|
||||
tmp<fvScalarMatrix> IiEq;
|
||||
|
||||
fvScalarMatrix IiEq
|
||||
if (!dom_.cacheDiv())
|
||||
{
|
||||
const surfaceScalarField Ji(dAve_ & mesh_.Sf());
|
||||
|
||||
IiEq =
|
||||
(
|
||||
fvm::div(Ji, ILambda_[lambdaI], "div(Ji,Ii_h)")
|
||||
+ fvm::Sp(k*omega_, ILambda_[lambdaI])
|
||||
==
|
||||
1.0/constant::mathematical::pi*omega_
|
||||
* (
|
||||
k*blackBody_.bLambda(lambdaI)
|
||||
+ absorptionEmission_.ECont(lambdaI)/4
|
||||
)
|
||||
);
|
||||
}
|
||||
else
|
||||
{
|
||||
IiEq =
|
||||
(
|
||||
dom_.fvRayDiv(myRayId_, lambdaI)
|
||||
+ fvm::Sp(k*omega_, ILambda_[lambdaI])
|
||||
==
|
||||
1.0/constant::mathematical::pi*omega_
|
||||
* (
|
||||
k*blackBody_.bLambda(lambdaI)
|
||||
+ absorptionEmission_.ECont(lambdaI)/4
|
||||
)
|
||||
);
|
||||
}
|
||||
|
||||
IiEq().relax();
|
||||
|
||||
const solverPerformance ILambdaSol = solve
|
||||
(
|
||||
fvm::div(Ji, ILambda_[lambdaI], "div(Ji,Ii_h)")
|
||||
+ fvm::Sp(k*omega_, ILambda_[lambdaI])
|
||||
==
|
||||
1.0/constant::mathematical::pi*omega_
|
||||
*(
|
||||
k*blackBody_.bLambda(lambdaI)
|
||||
+ absorptionEmission_.ECont(lambdaI)/4
|
||||
)
|
||||
IiEq(),
|
||||
mesh_.solver("Ii")
|
||||
);
|
||||
|
||||
IiEq.relax();
|
||||
const scalar initialRes =
|
||||
ILambdaSol.initialResidual()*omega_/dom_.omegaMax();
|
||||
|
||||
scalar eqnResidual = solve
|
||||
(
|
||||
IiEq,
|
||||
mesh_.solver("Ii")
|
||||
).initialResidual();
|
||||
|
||||
maxResidual = max(eqnResidual, maxResidual);
|
||||
maxResidual = max(initialRes, maxResidual);
|
||||
}
|
||||
|
||||
return maxResidual;
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -38,6 +38,7 @@ SourceFiles
|
||||
#include "absorptionEmissionModel.H"
|
||||
#include "blackBodyEmission.H"
|
||||
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
namespace Foam
|
||||
@ -111,6 +112,9 @@ private:
|
||||
//- Global ray id - incremented in constructor
|
||||
static label rayId;
|
||||
|
||||
//- My ray Id
|
||||
label myRayId_;
|
||||
|
||||
|
||||
// Private Member Functions
|
||||
|
||||
@ -209,6 +213,7 @@ public:
|
||||
|
||||
//- Return the radiative intensity for a given wavelength
|
||||
inline const volScalarField& ILambda(const label lambdaI) const;
|
||||
|
||||
};
|
||||
|
||||
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -252,8 +252,8 @@ void Foam::radiation::absorptionEmissionModel::correct
|
||||
PtrList<volScalarField>& aj
|
||||
) const
|
||||
{
|
||||
a.internalField() = this->a();
|
||||
aj[0].internalField() = a.internalField();
|
||||
a = this->a();
|
||||
aj[0] = a;
|
||||
}
|
||||
|
||||
|
||||
|
||||
Reference in New Issue
Block a user