49 void Foam::radiation::viewFactor::initialise()
51 const polyBoundaryMesh& coarsePatches = coarseMesh_.boundaryMesh();
52 const volScalarField::Boundary& Qrp = Qr_.boundaryField();
60 if ((isA<fixedValueFvPatchScalarField>(QrPatchi)))
62 selectedPatches_[count] = QrPatchi.patch().index();
63 nLocalCoarseFaces_ += coarsePatches[
patchi].size();
68 selectedPatches_.resize(count--);
72 Pout<<
"radiation::viewFactor::initialise() Selected patches:" 73 << selectedPatches_ <<
endl;
74 Pout<<
"radiation::viewFactor::initialise() Number of coarse faces:" 75 << nLocalCoarseFaces_ <<
endl;
78 totalNCoarseFaces_ = nLocalCoarseFaces_;
79 reduce(totalNCoarseFaces_, sumOp<label>());
84 <<
"Total number of clusters : " << totalNCoarseFaces_ <<
endl;
92 mesh_.facesInstance(),
105 mesh_.facesInstance(),
113 IOList<label> consMapDim
118 mesh_.facesInstance(),
131 Xfer<labelListList>(subMap),
132 Xfer<labelListList>(constructMap)
141 mesh_.facesInstance(),
154 mesh_.facesInstance(),
170 globalIndex globalNumbering(nLocalCoarseFaces_);
182 <<
"Insert elements in the matrix..." <<
endl;
191 globalFaceFacesProc[proci],
198 bool smoothing =
readBool(coeffs_.lookup(
"smoothing"));
204 <<
"Smoothing the matrix..." <<
endl;
207 for (
label i=0; i<totalNCoarseFaces_; i++)
210 for (
label j=0; j<totalNCoarseFaces_; j++)
212 sumF += Fmatrix_()(i, j);
214 scalar delta = sumF - 1.0;
215 for (
label j=0; j<totalNCoarseFaces_; j++)
217 Fmatrix_()(i, j) *= (1.0 - delta/(sumF + 0.001));
222 constEmissivity_ =
readBool(coeffs_.lookup(
"constantEmissivity"));
223 if (constEmissivity_)
230 pivotIndices_.setSize(CLU_().m());
246 mesh_.facesInstance(),
281 selectedPatches_(mesh_.
boundary().size(), -1),
282 totalNCoarseFaces_(0),
283 nLocalCoarseFaces_(0),
284 constEmissivity_(false),
292 Foam::radiation::viewFactor::viewFactor
317 mesh_.polyMesh::instance(),
340 totalNCoarseFaces_(0),
341 nLocalCoarseFaces_(0),
342 constEmissivity_(
false),
371 void Foam::radiation::viewFactor::insertMatrixElements
380 forAll(viewFactors, facei)
383 const labelList& globalFaces = globalFaceFaces[facei];
388 Fmatrix[globalI][globalFaces[i]] = vf[i];
399 scalarField compactCoarseT(map_->constructSize(), 0.0);
400 scalarField compactCoarseE(map_->constructSize(), 0.0);
401 scalarField compactCoarseHo(map_->constructSize(), 0.0);
412 forAll(selectedPatches_, i)
414 label patchID = selectedPatches_[i];
440 const labelList& agglom = finalAgglom_[patchID];
445 forAll(coarseToFine, coarseI)
447 const label coarseFaceID = coarsePatchFace[coarseI];
448 const labelList& fineFaces = coarseToFine[coarseFaceID];
454 scalar area =
sum(fineSf());
458 label facei = fineFaces[j];
459 Tave[coarseI] += (Tp[facei]*sf[facei])/area;
460 Eave[coarseI] += (eb[facei]*sf[facei])/area;
461 Hoiave[coarseI] += (Hoi[facei]*sf[facei])/area;
466 localCoarseTave.
append(Tave);
467 localCoarseEave.
append(Eave);
468 localCoarseHoave.
append(Hoiave);
474 SubList<scalar>(compactCoarseHo,nLocalCoarseFaces_) = localCoarseHoave;
477 map_->distribute(compactCoarseT);
478 map_->distribute(compactCoarseE);
479 map_->distribute(compactCoarseHo);
482 labelList compactGlobalIds(map_->constructSize(), 0.0);
484 labelList localGlobalIds(nLocalCoarseFaces_);
486 for(
label k = 0;
k < nLocalCoarseFaces_;
k++)
497 map_->distribute(compactGlobalIds);
507 T[compactGlobalIds[i]] = compactCoarseT[i];
508 E[compactGlobalIds[i]] = compactCoarseE[i];
509 QrExt[compactGlobalIds[i]] = compactCoarseHo[i];
526 if (!constEmissivity_)
530 for (
label i=0; i<totalNCoarseFaces_; i++)
532 for (
label j=0; j<totalNCoarseFaces_; j++)
534 scalar invEj = 1.0/E[j];
540 C(i, j) = invEj - (invEj - 1.0)*Fmatrix_()(i, j);
541 q[i] += (Fmatrix_()(i, j) - 1.0)*sigmaT4 - QrExt[j];
545 C(i, j) = (1.0 - invEj)*Fmatrix_()(i, j);
546 q[i] += Fmatrix_()(i, j)*sigmaT4;
552 Info<<
"\nSolving view factor equations..." <<
endl;
560 if (iterCounter_ == 0)
562 for (
label i=0; i<totalNCoarseFaces_; i++)
564 for (
label j=0; j<totalNCoarseFaces_; j++)
566 scalar invEj = 1.0/E[j];
569 CLU_()(i, j) = invEj-(invEj-1.0)*Fmatrix_()(i, j);
573 CLU_()(i, j) = (1.0 - invEj)*Fmatrix_()(i, j);
581 <<
"\nDecomposing C matrix..." <<
endl;
587 for (
label i=0; i<totalNCoarseFaces_; i++)
589 for (
label j=0; j<totalNCoarseFaces_; j++)
597 q[i] += (Fmatrix_()(i, j) - 1.0)*sigmaT4 - QrExt[j];
601 q[i] += Fmatrix_()(i, j)*sigmaT4;
609 <<
"\nLU Back substitute C matrix.." <<
endl;
621 label globCoarseId = 0;
622 forAll(selectedPatches_, i)
624 const label patchID = selectedPatches_[i];
630 const labelList& agglom = finalAgglom_[patchID];
638 scalar heatFlux = 0.0;
639 forAll(coarseToFine, coarseI)
643 const label coarseFaceID = coarsePatchFace[coarseI];
644 const labelList& fineFaces = coarseToFine[coarseFaceID];
647 label facei = fineFaces[
k];
649 Qrp[facei] = q[globalCoarse];
650 heatFlux += Qrp[facei]*sf[facei];
663 scalar heatFlux =
gSum(Qrp*magSf);
666 <<
"Total heat transfer rate at patch: "
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
label toGlobal(const label i) const
From local to global.
IOList< labelList > labelListIOList
Label container classes.
A list of keyword definitions, which are a keyword followed by any number of values (e...
const fileName & facesInstance() const
Return the current instance directory for faces.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
void LUDecompose(scalarSquareMatrix &matrix, labelList &pivotIndices)
LU decompose the matrix with pivoting.
labelListList invertOneToMany(const label len, const labelUList &)
Invert one-to-many map. Unmapped elements will be size 0.
To & refCast(From &r)
Reference type cast template function.
const volScalarField & T_
Reference to the temperature field.
void size(const label)
Override size to be inconsistent with allocated storage.
IOobject(const word &name, const fileName &instance, const objectRegistry ®istry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
static bool master(const label communicator=0)
Am I the master process.
scalarField emissivity() const
Calculate corresponding emissivity field.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
const Type & value() const
Return const reference to value.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
label k
Boltzmann constant.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Macros for easy insertion into run-time selection tables.
const labelListList & patchFaceMap() const
From patchFace on this back to original mesh or agglomeration.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
A List obtained as a section of another List.
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
Type gSum(const FieldField< Field, Type > &f)
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
virtual bool read()=0
Read radiationProperties dictionary.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
fvPatchField< scalar > fvPatchScalarField
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
const fvMesh & mesh_
Reference to the mesh database.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Top level model for radiation modelling.
faceListList boundary(nPatches)
void calculate()
Solve system of equation(s)
prefixOSstream Pout(cout,"Pout")
const scalarList & Qro()
Return external radiative heat flux.
virtual ~viewFactor()
Destructor.
defineTypeNameAndDebug(combustionModel, 0)
bool read()
Read radiation properties dictionary.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
volScalarField sf(fieldObject, mesh)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
virtual tmp< DimensionedField< scalar, volMesh > > Ru() const
Source term component (constant)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
This boundary condition provides a grey-diffuse condition for radiative heat flux, Qr, for use with the view factor model.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
const dimensionedScalar F
Faraday constant: default SI units: [C/mol].
void relax(const scalar alpha)
Relax field (for steady-state solution).
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A List with indirect addressing.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
IOList< scalarList > scalarListIOList
Scalar container classes.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Mesh consisting of general polyhedral cells.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
A class for managing temporary objects.
A patch is a list of labels that address the faces in the global face list.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
#define addToRadiationRunTimeSelectionTables(model)
SquareMatrix< scalar > scalarSquareMatrix
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
void storePrevIter() const
Store the field as the previous iteration value.
const word & name() const
Return reference to name.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
const Time & time() const
Return the top-level database.
void LUsolve(scalarSquareMatrix &matrix, List< Type > &source)
Solve the matrix using LU decomposition with pivoting.
void LUBacksubstitute(const scalarSquareMatrix &luMmatrix, const labelList &pivotIndices, List< Type > &source)
LU back-substitution with given source, returning the solution.
label size() const
Return the number of elements in the UPtrList.
#define InfoInFunction
Report an information message using Foam::Info.