Go to the documentation of this file.
59 selectedPatches_ = mesh_.boundaryMesh().
indices(viewFactorWalls);
61 for (
const label patchi : selectedPatches_)
63 nLocalCoarseFaces_ += coarsePatches[patchi].size();
68 Pout<<
"radiation::viewFactor::initialise() Selected patches:"
69 << selectedPatches_ <<
endl;
70 Pout<<
"radiation::viewFactor::initialise() Number of coarse faces:"
71 << nLocalCoarseFaces_ <<
endl;
74 totalNCoarseFaces_ = nLocalCoarseFaces_;
80 <<
"Total number of clusters : " << totalNCoarseFaces_ <<
endl;
90 mesh_.facesInstance(),
104 mesh_.facesInstance(),
117 mesh_.facesInstance(),
145 <<
"Insert elements in the matrix..." <<
endl;
154 globalFaceFacesProc[procI],
161 if (coeffs_.get<
bool>(
"smoothing"))
166 <<
"Smoothing the matrix..." <<
endl;
169 for (
label i=0; i<totalNCoarseFaces_; i++)
172 for (
label j=0; j<totalNCoarseFaces_; j++)
174 sumF += Fmatrix_()(i, j);
177 const scalar
delta = sumF - 1.0;
178 for (
label j=0; j<totalNCoarseFaces_; j++)
180 Fmatrix_()(i, j) *= (1.0 -
delta/(sumF + 0.001));
185 coeffs_.readEntry(
"constantEmissivity", constEmissivity_);
186 if (constEmissivity_)
193 pivotIndices_.setSize(CLU_().m());
197 coeffs_.readIfPresent(
"useSolarLoad", useSolarLoad_);
201 const dictionary& solarDict = this->subDict(
"solarLoadCoeffs");
202 solarLoad_.reset(
new solarLoad(solarDict, T_));
204 if (solarLoad_->nBands() != nBands_)
207 <<
"Solar radiation and view factor band numbers "
212 Info<<
"Creating Solar Load Model " <<
nl;
227 mesh_.facesInstance(),
239 "coarse:" + mesh_.
name(),
263 selectedPatches_(mesh_.
boundary().size(), -1),
264 totalNCoarseFaces_(0),
265 nLocalCoarseFaces_(0),
266 constEmissivity_(false),
269 useSolarLoad_(false),
271 nBands_(coeffs_.lookupOrDefault<
label>(
"nBands", 1))
277 Foam::radiation::viewFactor::viewFactor
289 mesh_.facesInstance(),
301 "coarse:" + mesh_.name(),
302 mesh_.polyMesh::instance(),
316 mesh_.time().timeName(),
325 selectedPatches_(mesh_.boundary().size(), -1),
326 totalNCoarseFaces_(0),
327 nLocalCoarseFaces_(0),
328 constEmissivity_(
false),
331 useSolarLoad_(
false),
333 nBands_(coeffs_.lookupOrDefault<
label>(
"nBands", 1))
361 forAll(viewFactors, faceI)
364 const labelList& globalFaces = globalFaceFaces[faceI];
369 Fmatrix[globalI][globalFaces[i]] = vf[i];
382 solarLoad_->calculate();
387 volScalarField::Boundary& qrBf = qr_.boundaryFieldRef();
394 for (
label bandI = 0; bandI < nBands_; bandI++)
396 scalarField compactCoarseT4(map_->constructSize(), 0.0);
397 scalarField compactCoarseE(map_->constructSize(), 0.0);
398 scalarField compactCoarseHo(map_->constructSize(), 0.0);
405 forAll(selectedPatches_, i)
431 coarseMesh_.patchFaceMap()[
patchID];
444 forAll(coarseToFine, coarseI)
446 const label coarseFaceID = coarsePatchFace[coarseI];
447 const labelList& fineFaces = coarseToFine[coarseFaceID];
454 const scalar
area =
sum(fineSf());
459 label facei = fineFaces[j];
460 T4ave[coarseI] += (
pow4(Tp[facei])*sf[facei])/
area;
461 Eave[coarseI] += (eb[facei]*sf[facei])/
area;
462 Hoiave[coarseI] += (Hoi[facei]*sf[facei])/
area;
467 localCoarseT4ave.
append(T4ave);
468 localCoarseEave.
append(Eave);
469 localCoarseHoave.
append(Hoiave);
481 map_->distribute(compactCoarseT4);
482 map_->distribute(compactCoarseE);
483 map_->distribute(compactCoarseHo);
498 map_->distribute(compactGlobalIds);
506 forAll(compactCoarseT4, i)
508 T4[compactGlobalIds[i]] = compactCoarseT4[i];
509 E[compactGlobalIds[i]] = compactCoarseE[i];
510 qrExt[compactGlobalIds[i]] = compactCoarseHo[i];
524 if (!constEmissivity_)
528 for (
label i=0; i<totalNCoarseFaces_; i++)
530 for (
label j=0; j<totalNCoarseFaces_; j++)
532 const scalar invEj = 1.0/E[j];
533 const scalar sigmaT4 =
538 C(i, j) = invEj - (invEj - 1.0)*Fmatrix_()(i, j);
540 (Fmatrix_()(i, j) - 1.0)*sigmaT4 + qrExt[j];
544 C(i, j) = (1.0 - invEj)*Fmatrix_()(i, j);
545 q[i] += Fmatrix_()(i, j)*sigmaT4;
551 Info<<
"Solving view factor equations for band :"
560 if (iterCounter_ == 0)
562 for (
label i=0; i<totalNCoarseFaces_; i++)
564 for (
label j=0; j<totalNCoarseFaces_; j++)
566 const scalar invEj = 1.0/E[j];
570 invEj-(invEj-1.0)*Fmatrix_()(i, j);
574 CLU_()(i, j) = (1.0 - invEj)*Fmatrix_()(i, j);
582 <<
"\nDecomposing C matrix..." <<
endl;
588 for (
label i=0; i<totalNCoarseFaces_; i++)
590 for (
label j=0; j<totalNCoarseFaces_; j++)
592 const scalar sigmaT4 =
598 (Fmatrix_()(i, j) - 1.0)*sigmaT4 + qrExt[j];
602 q[i] += Fmatrix_()(i, j)*sigmaT4;
608 Info<<
"Solving view factor equations for band : "
622 label globCoarseId = 0;
637 coarseMesh_.patchFaceMap()[
patchID];
639 scalar heatFlux = 0.0;
640 forAll(coarseToFine, coarseI)
646 const label coarseFaceID = coarsePatchFace[coarseI];
647 const labelList& fineFaces = coarseToFine[coarseFaceID];
650 label faceI = fineFaces[
k];
652 qrp[faceI] = q[globalCoarse];
653 heatFlux += qrp[faceI]*sf[faceI];
666 const scalar heatFlux =
gSum(qrp*magSf);
669 <<
"Total heat transfer rate at patch: "
687 mesh_.time().timeName(),
710 mesh_.time().timeName(),
int debug
Static debugging option.
void LUsolve(scalarSquareMatrix &matrix, List< Type > &source)
Defines the attributes of an object for which implicit objectRegistry management is supported,...
#define InfoInFunction
Report an information message using Foam::Info.
A class for handling words, derived from Foam::string.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
label localStart() const
My local start.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero.
Different types of constants.
IOmapDistribute is derived from mapDistribute and IOobject to give the mapDistribute automatic IO fun...
A List obtained as a section of another List.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
static const boundaryRadiationProperties & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
volVectorField F(fluid.F())
label localSize() const
My local size.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
Type gSum(const FieldField< Field, Type > &f)
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
void insertMatrixElements(const globalIndex &index, const label fromProci, const labelListList &globalFaceFaces, const scalarListList &viewFactors, scalarSquareMatrix &matrix)
Insert view factors into main matrix.
Mesh consisting of general polyhedral cells.
tmp< scalarField > qro(label bandI=0) const
Return external + solar load radiative heat flux.
#define forAll(list, i)
Loop across all elements in list.
tmp< scalarField > emissivity(const label patchI, const label bandI=0, vectorField *incomingDirection=nullptr, scalarField *T=nullptr) const
Access boundary emissivity on patch.
void initialise()
Initialise.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
dimensionedScalar pow4(const dimensionedScalar &ds)
virtual bool read()=0
Read radiationProperties dictionary.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
dimensionedScalar pow3(const dimensionedScalar &ds)
messageStream Info
Information stream (uses stdout - output is on the master only)
A patch is a list of labels that address the faces in the global face list.
word name(const complex &c)
Return string representation of complex.
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
Boundary radiation properties holder.
SquareMatrix< scalar > scalarSquareMatrix
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
void calculate()
Solve system of equation(s)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
errorManip< error > abort(error &err)
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
To & refCast(From &r)
Reference type cast template function.
void LUBacksubstitute(const scalarSquareMatrix &luMmatrix, const labelList &pivotIndices, List< Type > &source)
labelList indices(const keyType &key, const bool useGroups=true) const
Return patch indices for all matches.
void LUDecompose(scalarSquareMatrix &matrix, labelList &pivotIndices)
LU decompose the matrix with pivoting.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
static bool master(const label communicator=0)
Am I the master process.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
The solar load radiation model includes Sun primary hits, their reflective fluxes and diffusive sky r...
This boundary condition provides a grey-diffuse condition for radiative heat flux,...
Top level model for radiation modelling.
label k
Boltzmann constant.
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
A List of objects of type <T> with automated input and output.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
static const word viewFactorWalls
Static name for view factor walls.
const wordList area
Standard area field types (scalar, vector, tensor, etc)
A List with indirect addressing.
virtual tmp< volScalarField::Internal > Ru() const
Source term component (constant)
bool read()
Read radiation properties dictionary.
#define addToRadiationRunTimeSelectionTables(model)
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
Graphite solid properties.
labelListList invertOneToMany(const label len, const labelUList &map)
Invert one-to-many map. Unmapped elements will be size 0.
defineTypeNameAndDebug(combustionModel, 0)
label toGlobal(const label i) const
From local to global index.