viewFactor.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "viewFactor.H"
27 #include "surfaceFields.H"
28 #include "constants.H"
30 #include "typeInfo.H"
32 
33 using namespace Foam::constant;
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace radiationModels
40 {
41  defineTypeNameAndDebug(viewFactor, 0);
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
48 
49 void Foam::radiationModels::viewFactor::initialise()
50 {
51  const polyBoundaryMesh& coarsePatches = coarseMesh_.boundaryMesh();
52  const volScalarField::Boundary& qrp = qr_.boundaryField();
53 
54  label count = 0;
55  forAll(qrp, patchi)
56  {
57  // const polyPatch& pp = mesh_.boundaryMesh()[patchi];
58  const fvPatchScalarField& qrPatchi = qrp[patchi];
59 
60  if ((isA<fixedValueFvPatchScalarField>(qrPatchi)))
61  {
62  selectedPatches_[count] = qrPatchi.patch().index();
63  nLocalCoarseFaces_ += coarsePatches[patchi].size();
64  count++;
65  }
66  }
67 
68  selectedPatches_.resize(count--);
69 
70  if (debug)
71  {
72  Pout<< "radiationModels::viewFactor::initialise() "
73  << "Selected patches:" << selectedPatches_ << endl;
74  Pout<< "radiationModels::viewFactor::initialise() "
75  << "Number of coarse faces:" << nLocalCoarseFaces_ << endl;
76  }
77 
78  totalNCoarseFaces_ = nLocalCoarseFaces_;
79  reduce(totalNCoarseFaces_, sumOp<label>());
80 
81  if (debug && Pstream::master())
82  {
84  << "Total number of clusters : " << totalNCoarseFaces_ << endl;
85  }
86 
87  labelListIOList subMap
88  (
89  IOobject
90  (
91  "subMap",
92  mesh_.facesInstance(),
93  mesh_,
96  false
97  )
98  );
99 
100  labelListIOList constructMap
101  (
102  IOobject
103  (
104  "constructMap",
105  mesh_.facesInstance(),
106  mesh_,
109  false
110  )
111  );
112 
113  IOList<label> consMapDim
114  (
115  IOobject
116  (
117  "constructMapDim",
118  mesh_.facesInstance(),
119  mesh_,
122  false
123  )
124  );
125 
126  map_.reset
127  (
128  new distributionMap
129  (
130  consMapDim[0],
131  move(subMap),
132  move(constructMap)
133  )
134  );
135 
136  scalarListIOList FmyProc
137  (
138  IOobject
139  (
140  "F",
141  mesh_.facesInstance(),
142  mesh_,
145  false
146  )
147  );
148 
149  labelListIOList globalFaceFaces
150  (
151  IOobject
152  (
153  "globalFaceFaces",
154  mesh_.facesInstance(),
155  mesh_,
158  false
159  )
160  );
161 
162  List<labelListList> globalFaceFacesProc(Pstream::nProcs());
163  globalFaceFacesProc[Pstream::myProcNo()] = globalFaceFaces;
164  Pstream::gatherList(globalFaceFacesProc);
165 
166  List<scalarListList> F(Pstream::nProcs());
167  F[Pstream::myProcNo()] = FmyProc;
169 
170  globalIndex globalNumbering(nLocalCoarseFaces_);
171 
172  if (Pstream::master())
173  {
174  Fmatrix_.reset
175  (
176  new scalarSquareMatrix(totalNCoarseFaces_, 0.0)
177  );
178 
179  if (debug)
180  {
182  << "Insert elements in the matrix..." << endl;
183  }
184 
185  for (label proci = 0; proci < Pstream::nProcs(); proci++)
186  {
187  insertMatrixElements
188  (
189  globalNumbering,
190  proci,
191  globalFaceFacesProc[proci],
192  F[proci],
193  Fmatrix_()
194  );
195  }
196 
197 
198  bool smoothing = readBool(coeffs_.lookup("smoothing"));
199  if (smoothing)
200  {
201  if (debug)
202  {
204  << "Smoothing the matrix..." << endl;
205  }
206 
207  for (label i=0; i<totalNCoarseFaces_; i++)
208  {
209  scalar sumF = 0.0;
210  for (label j=0; j<totalNCoarseFaces_; j++)
211  {
212  sumF += Fmatrix_()(i, j);
213  }
214 
215  const scalar delta = sumF - 1.0;
216  for (label j=0; j<totalNCoarseFaces_; j++)
217  {
218  Fmatrix_()(i, j) *= (1.0 - delta/(sumF + 0.001));
219  }
220  }
221  }
222 
223  constEmissivity_ = readBool(coeffs_.lookup("constantEmissivity"));
224  if (constEmissivity_)
225  {
226  CLU_.reset
227  (
228  new scalarSquareMatrix(totalNCoarseFaces_, 0.0)
229  );
230 
231  pivotIndices_.setSize(CLU_().m());
232  }
233  }
234 }
235 
236 
237 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
238 
240 :
241  radiationModel(typeName, T),
242  finalAgglom_
243  (
244  IOobject
245  (
246  "finalAgglom",
247  mesh_.facesInstance(),
248  mesh_,
249  IOobject::MUST_READ,
250  IOobject::NO_WRITE,
251  false
252  )
253  ),
254  map_(),
255  coarseMesh_
256  (
257  IOobject
258  (
259  mesh_.name(),
260  mesh_.polyMesh::instance(),
261  mesh_.time(),
262  IOobject::NO_READ,
263  IOobject::NO_WRITE
264  ),
265  mesh_,
266  finalAgglom_
267  ),
268  qr_
269  (
270  IOobject
271  (
272  "qr",
273  mesh_.time().timeName(),
274  mesh_,
275  IOobject::MUST_READ,
276  IOobject::AUTO_WRITE
277  ),
278  mesh_
279  ),
280  Fmatrix_(),
281  CLU_(),
282  selectedPatches_(mesh_.boundary().size(), -1),
283  totalNCoarseFaces_(0),
284  nLocalCoarseFaces_(0),
285  constEmissivity_(false),
286  iterCounter_(0),
287  pivotIndices_(0)
288 {
289  initialise();
290 }
291 
292 
294 (
295  const dictionary& dict,
296  const volScalarField& T
297 )
298 :
299  radiationModel(typeName, dict, T),
300  finalAgglom_
301  (
302  IOobject
303  (
304  "finalAgglom",
306  mesh_,
309  false
310  )
311  ),
312  map_(),
313  coarseMesh_
314  (
315  IOobject
316  (
317  mesh_.name(),
318  mesh_.polyMesh::instance(),
319  mesh_.time(),
322  ),
323  mesh_,
324  finalAgglom_
325  ),
326  qr_
327  (
328  IOobject
329  (
330  "qr",
331  mesh_.time().timeName(),
332  mesh_,
335  ),
336  mesh_
337  ),
338  Fmatrix_(),
339  CLU_(),
340  selectedPatches_(mesh_.boundary().size(), -1),
341  totalNCoarseFaces_(0),
342  nLocalCoarseFaces_(0),
343  constEmissivity_(false),
344  iterCounter_(0),
345  pivotIndices_(0)
346 {
347  initialise();
348 }
349 
350 
351 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
352 
354 {}
355 
356 
357 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
358 
360 {
361  if (radiationModel::read())
362  {
363  return true;
364  }
365  else
366  {
367  return false;
368  }
369 }
370 
371 
372 void Foam::radiationModels::viewFactor::insertMatrixElements
373 (
374  const globalIndex& globalNumbering,
375  const label proci,
376  const labelListList& globalFaceFaces,
377  const scalarListList& viewFactors,
378  scalarSquareMatrix& Fmatrix
379 )
380 {
381  forAll(viewFactors, facei)
382  {
383  const scalarList& vf = viewFactors[facei];
384  const labelList& globalFaces = globalFaceFaces[facei];
385 
386  label globalI = globalNumbering.toGlobal(proci, facei);
387  forAll(globalFaces, i)
388  {
389  Fmatrix[globalI][globalFaces[i]] = vf[i];
390  }
391  }
392 }
393 
394 
396 {
397  // Store previous iteration
398  qr_.storePrevIter();
399 
400  scalarField compactCoarseT4(map_->constructSize(), 0.0);
401  scalarField compactCoarseE(map_->constructSize(), 0.0);
402  scalarField compactCoarseHo(map_->constructSize(), 0.0);
403 
404  globalIndex globalNumbering(nLocalCoarseFaces_);
405 
406  // Fill local averaged(T), emissivity(E) and external heatFlux(Ho)
407  DynamicList<scalar> localCoarseT4ave(nLocalCoarseFaces_);
408  DynamicList<scalar> localCoarseEave(nLocalCoarseFaces_);
409  DynamicList<scalar> localCoarseHoave(nLocalCoarseFaces_);
410 
412 
413  forAll(selectedPatches_, i)
414  {
415  label patchID = selectedPatches_[i];
416 
417  const scalarField& Tp = T_.boundaryField()[patchID];
418  const scalarField& sf = mesh_.magSf().boundaryField()[patchID];
419 
420  fvPatchScalarField& qrPatch = qrBf[patchID];
421 
423  refCast
424  <
426  >(qrPatch);
427 
428  const scalarList eb = qrp.emissivity();
429 
430  const scalarList& Hoi = qrp.qro();
431 
432  const polyPatch& pp = coarseMesh_.boundaryMesh()[patchID];
433  const labelList& coarsePatchFace = coarseMesh_.patchFaceMap()[patchID];
434 
435  scalarList T4ave(pp.size(), 0.0);
436  scalarList Eave(pp.size(), 0.0);
437  scalarList Hoiave(pp.size(), 0.0);
438 
439  if (pp.size() > 0)
440  {
441  const labelList& agglom = finalAgglom_[patchID];
442  label nAgglom = max(agglom) + 1;
443 
444  labelListList coarseToFine(invertOneToMany(nAgglom, agglom));
445 
446  forAll(coarseToFine, coarseI)
447  {
448  const label coarseFaceID = coarsePatchFace[coarseI];
449  const labelList& fineFaces = coarseToFine[coarseFaceID];
450  UIndirectList<scalar> fineSf
451  (
452  sf,
453  fineFaces
454  );
455 
456  const scalar area = sum(fineSf());
457 
458  // Temperature, emissivity and external flux area weighting
459  forAll(fineFaces, j)
460  {
461  label facei = fineFaces[j];
462  T4ave[coarseI] += (pow4(Tp[facei])*sf[facei])/area;
463  Eave[coarseI] += (eb[facei]*sf[facei])/area;
464  Hoiave[coarseI] += (Hoi[facei]*sf[facei])/area;
465  }
466  }
467  }
468 
469  localCoarseT4ave.append(T4ave);
470  localCoarseEave.append(Eave);
471  localCoarseHoave.append(Hoiave);
472  }
473 
474  // Fill the local values to distribute
475  SubList<scalar>(compactCoarseT4, nLocalCoarseFaces_) = localCoarseT4ave;
476  SubList<scalar>(compactCoarseE, nLocalCoarseFaces_) = localCoarseEave;
477  SubList<scalar>(compactCoarseHo, nLocalCoarseFaces_) = localCoarseHoave;
478 
479  // Distribute data
480  map_->distribute(compactCoarseT4);
481  map_->distribute(compactCoarseE);
482  map_->distribute(compactCoarseHo);
483 
484  // Distribute local global ID
485  labelList compactGlobalIds(map_->constructSize(), 0.0);
486 
487  labelList localGlobalIds(nLocalCoarseFaces_);
488 
489  for(label k = 0; k < nLocalCoarseFaces_; k++)
490  {
491  localGlobalIds[k] = globalNumbering.toGlobal(Pstream::myProcNo(), k);
492  }
493 
495  (
496  compactGlobalIds,
497  nLocalCoarseFaces_
498  ) = localGlobalIds;
499 
500  map_->distribute(compactGlobalIds);
501 
502  // Create global size vectors
503  scalarField T4(totalNCoarseFaces_, 0.0);
504  scalarField E(totalNCoarseFaces_, 0.0);
505  scalarField qrExt(totalNCoarseFaces_, 0.0);
506 
507  // Fill lists from compact to global indexes.
508  forAll(compactCoarseT4, i)
509  {
510  T4[compactGlobalIds[i]] = compactCoarseT4[i];
511  E[compactGlobalIds[i]] = compactCoarseE[i];
512  qrExt[compactGlobalIds[i]] = compactCoarseHo[i];
513  }
514 
518 
522 
523  // Net radiation
524  scalarField q(totalNCoarseFaces_, 0.0);
525 
526  if (Pstream::master())
527  {
528  // Variable emissivity
529  if (!constEmissivity_)
530  {
531  scalarSquareMatrix C(totalNCoarseFaces_, 0.0);
532 
533  for (label i=0; i<totalNCoarseFaces_; i++)
534  {
535  for (label j=0; j<totalNCoarseFaces_; j++)
536  {
537  const scalar invEj = 1.0/E[j];
538  const scalar sigmaT4 = physicoChemical::sigma.value()*T4[j];
539 
540  if (i==j)
541  {
542  C(i, j) = invEj - (invEj - 1.0)*Fmatrix_()(i, j);
543  q[i] += (Fmatrix_()(i, j) - 1.0)*sigmaT4 - qrExt[j];
544  }
545  else
546  {
547  C(i, j) = (1.0 - invEj)*Fmatrix_()(i, j);
548  q[i] += Fmatrix_()(i, j)*sigmaT4;
549  }
550 
551  }
552  }
553 
554  Info<< "\nSolving view factor equations..." << endl;
555 
556  // Negative coming into the fluid
557  LUsolve(C, q);
558  }
559  else // Constant emissivity
560  {
561  // Initial iter calculates CLU and chaches it
562  if (iterCounter_ == 0)
563  {
564  for (label i=0; i<totalNCoarseFaces_; i++)
565  {
566  for (label j=0; j<totalNCoarseFaces_; j++)
567  {
568  const scalar invEj = 1.0/E[j];
569  if (i==j)
570  {
571  CLU_()(i, j) = invEj-(invEj-1.0)*Fmatrix_()(i, j);
572  }
573  else
574  {
575  CLU_()(i, j) = (1.0 - invEj)*Fmatrix_()(i, j);
576  }
577  }
578  }
579 
580  if (debug)
581  {
583  << "\nDecomposing C matrix..." << endl;
584  }
585 
586  LUDecompose(CLU_(), pivotIndices_);
587  }
588 
589  for (label i=0; i<totalNCoarseFaces_; i++)
590  {
591  for (label j=0; j<totalNCoarseFaces_; j++)
592  {
593  const scalar sigmaT4 =
595 
596  if (i==j)
597  {
598  q[i] += (Fmatrix_()(i, j) - 1.0)*sigmaT4 - qrExt[j];
599  }
600  else
601  {
602  q[i] += Fmatrix_()(i, j)*sigmaT4;
603  }
604  }
605  }
606 
607  if (debug)
608  {
610  << "\nLU Back substitute C matrix.." << endl;
611  }
612 
613  LUBacksubstitute(CLU_(), pivotIndices_, q);
614  iterCounter_ ++;
615  }
616  }
617 
618  // Scatter q and fill qr
621 
622  label globCoarseId = 0;
623  forAll(selectedPatches_, i)
624  {
625  const label patchID = selectedPatches_[i];
626  const polyPatch& pp = mesh_.boundaryMesh()[patchID];
627  if (pp.size() > 0)
628  {
629  scalarField& qrp = qrBf[patchID];
630  const scalarField& sf = mesh_.magSf().boundaryField()[patchID];
631  const labelList& agglom = finalAgglom_[patchID];
632  label nAgglom = max(agglom)+1;
633 
634  labelListList coarseToFine(invertOneToMany(nAgglom, agglom));
635 
636  const labelList& coarsePatchFace =
637  coarseMesh_.patchFaceMap()[patchID];
638 
639  scalar heatFlux = 0.0;
640  forAll(coarseToFine, coarseI)
641  {
642  label globalCoarse =
643  globalNumbering.toGlobal(Pstream::myProcNo(), globCoarseId);
644  const label coarseFaceID = coarsePatchFace[coarseI];
645  const labelList& fineFaces = coarseToFine[coarseFaceID];
646  forAll(fineFaces, k)
647  {
648  label facei = fineFaces[k];
649 
650  qrp[facei] = q[globalCoarse];
651  heatFlux += qrp[facei]*sf[facei];
652  }
653  globCoarseId ++;
654  }
655  }
656  }
657 
658  if (debug)
659  {
660  forAll(qrBf, patchID)
661  {
662  const scalarField& qrp = qrBf[patchID];
663  const scalarField& magSf = mesh_.magSf().boundaryField()[patchID];
664  const scalar heatFlux = gSum(qrp*magSf);
665 
667  << "Total heat transfer rate at patch: "
668  << patchID << " "
669  << heatFlux << endl;
670  }
671  }
672 
673  // Relax qr if necessary
674  qr_.relax();
675 }
676 
677 
679 {
680  return volScalarField::New
681  (
682  "Rp",
683  mesh_,
685  (
687  0
688  )
689  );
690 }
691 
692 
695 {
697  (
698  "Ru",
699  mesh_,
701  );
702 }
703 
704 // ************************************************************************* //
Collection of constants.
Foam::surfaceFields.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:453
This boundary condition provides a grey-diffuse condition for radiative heat flux, qr, for use with the view factor model.
Graphite solid properties.
Definition: C.H:48
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.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &)
Return a temporary field constructed from name, mesh.
const labelListList & patchFaceMap() const
From patchFace on this back to original mesh or agglomeration.
const dimensionedScalar F
Faraday constant: default SI units: [C/mol].
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:888
IOList< labelList > labelListIOList
Label container classes.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const Boundary & boundaryField() const
Return const-reference to the boundary field.
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.
Definition: ListOps.C:67
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:106
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
label count(const ListType &l, typename ListType::const_reference x)
Count the number of occurrences of a value in a list.
virtual bool read()=0
Read radiationProperties dictionary.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
scalarField emissivity() const
Calculate corresponding emissivity field.
label k
Boltzmann constant.
bool readBool(Istream &)
Definition: boolIO.C:60
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
void calculate()
Solve system of equation(s)
Definition: viewFactor.C:395
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Macros for easy insertion into run-time selection tables.
const dimensionSet dimLength
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const volScalarField & T_
Reference to the temperature field.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
const dimensionSet dimTime
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
A List obtained as a section of another List.
Definition: SubList.H:53
Type gSum(const FieldField< Field, Type > &f)
fvPatchField< scalar > fvPatchScalarField
static word timeName(const scalar, const int precision=curPrecision_)
Return time name of given scalar time.
Definition: Time.C:666
Top level model for radiation modelling.
virtual tmp< volScalarField::Internal > Ru() const
Source term component (constant)
Definition: viewFactor.C:694
virtual ~viewFactor()
Destructor.
Definition: viewFactor.C:353
bool read()
Read radiation properties dictionary.
Definition: viewFactor.C:359
const Type & value() const
Return const reference to value.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
word timeName
Definition: getTimeIndex.H:3
faceListList boundary(nPatches)
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
defineTypeNameAndDebug(combustionModel, 0)
const dimensionSet dimMass
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.
Definition: complex.C:47
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
viewFactor(const volScalarField &T)
Construct from components.
Definition: viewFactor.C:239
dimensionedScalar pow3(const dimensionedScalar &ds)
label toGlobal(const label i) const
From local to global.
Definition: globalIndexI.H:82
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
label patchi
void relax(const scalar alpha)
Relax field (for steady-state solution).
const word & name() const
Return reference to name.
Definition: fvMesh.H:386
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
radiationModel(const volScalarField &T)
Null constructor.
dimensionedScalar pow4(const dimensionedScalar &ds)
IOList< scalarList > scalarListIOList
Scalar container classes.
messageStream Info
const fvMesh & mesh_
Reference to the mesh database.
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
Definition: viewFactor.C:678
void storePrevIter() const
Store the field as the previous iteration value.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
A class for managing temporary objects.
Definition: PtrList.H:53
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
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 dimensionSet dimTemperature
Namespace for OpenFOAM.
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.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:800
#define InfoInFunction
Report an information message using Foam::Info.