laserDTRM.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017-2022 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "laserDTRM.H"
29#include "fvmLaplacian.H"
30#include "fvmSup.H"
32#include "scatterModel.H"
33#include "constants.H"
34#include "unitConversion.H"
35#include "interpolationCell.H"
37#include "Random.H"
38#include "OBJstream.H"
40
41using namespace Foam::constant;
42
43// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
44
45namespace Foam
46{
47 namespace radiation
48 {
51 }
52
54 (
56 "DTRMCloud",
57 0
58 );
59
60} // End namespace Foam
61
62
64Foam::radiation::laserDTRM::powerDistNames_
65{
66 { powerDistributionMode::pdGaussian, "Gaussian" },
67 { powerDistributionMode::pdManual, "manual" },
68 { powerDistributionMode::pdUniform, "uniform" },
69 { powerDistributionMode::pdGaussianPeak, "GaussianPeak" },
70};
71
72
73// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
74
75Foam::scalar Foam::radiation::laserDTRM::calculateIp(scalar r, scalar theta)
76{
77 const scalar t = mesh_.time().value();
78 const scalar power = laserPower_->value(t);
79 switch (mode_)
80 {
81 case pdGaussianPeak:
82 {
83 return I0_*exp(-2.0*sqr(r)/sqr(sigma_));
84 break;
85 }
86 case pdGaussian:
87 {
88 scalar I0 = power/(mathematical::twoPi*sqr(sigma_));
89
90 return I0*exp(-sqr(r)/2.0/sqr(sigma_));
91 break;
92 }
93 case pdManual:
94 {
95 return power*powerDistribution_()(theta, r);
96 break;
97 }
98 case pdUniform:
99 {
100 return power/(mathematical::pi*sqr(focalLaserRadius_));
101 break;
102 }
103 default:
104 {
106 << "Unhandled type " << powerDistNames_[mode_]
107 << abort(FatalError);
108 break;
109 }
110 }
111
112 return 0;
113}
114
115
116Foam::tmp<Foam::volVectorField> Foam::radiation::laserDTRM::nHatfv
117(
118 const volScalarField& alpha1,
120) const
121{
122 const dimensionedScalar deltaN
123 (
124 "deltaN",
125 1e-7/cbrt(average(mesh_.V()))
126 );
127
128 const volVectorField gradAlphaf
129 (
132 );
133
134 // Face unit interface normal
135 return gradAlphaf/(mag(gradAlphaf)+ deltaN);
136}
137
138
139void Foam::radiation::laserDTRM::initialiseReflection()
140{
141 if (found("reflectionModel"))
142 {
143 dictTable modelDicts(lookup("reflectionModel"));
144
145 forAllConstIters(modelDicts, iter)
146 {
147 const phasePairKey& key = iter.key();
148
149 reflections_.insert
150 (
151 key,
153 (
154 iter.val(),
155 mesh_
156 )
157 );
158 }
159
160 if (reflections_.size())
161 {
162 reflectionSwitch_ = true;
163 }
164
165 reflectionSwitch_ = returnReduce(reflectionSwitch_, orOp<bool>());
166 }
167}
168
169
170void Foam::radiation::laserDTRM::initialise()
171{
172 // Initialise the DTRM particles
173 DTRMCloud_.clear();
174
175 const scalar t = mesh_.time().value();
176 const vector lPosition = focalLaserPosition_->value(t);
177 const vector lDir = normalised(laserDirection_->value(t));
178
180 << "Laser position : " << lPosition << nl
181 << "Laser direction : " << lDir << endl;
182
183 // Find a vector on the area plane. Normal to laser direction
184 vector rArea = Zero;
185 scalar magr = 0.0;
186
187 {
188 Random rnd(1234);
189
190 while (magr < VSMALL)
191 {
192 vector v = rnd.sample01<vector>();
193 rArea = v - (v & lDir)*lDir;
194 magr = mag(rArea);
195 }
196 }
197 rArea.normalise();
198
199 scalar dr = focalLaserRadius_/ndr_;
200 scalar dTheta = mathematical::twoPi/ndTheta_;
201
202 nParticles_ = ndr_*ndTheta_;
203
204 switch (mode_)
205 {
206 case pdGaussianPeak:
207 {
208 I0_ = get<scalar>("I0");
209 sigma_ = get<scalar>("sigma");
210 break;
211 }
212 case pdGaussian:
213 {
214 sigma_ = get<scalar>("sigma");
215 break;
216 }
217 case pdManual:
218 {
219 powerDistribution_.reset
220 (
222 );
223
224 break;
225 }
226 case pdUniform:
227 {
228 break;
229 }
230 }
231
232 // Count the number of missed positions
233 label nMissed = 0;
234
235 // Target position
236 point p1 = vector::zero;
237
238 // Seed DTRM particles
239 // TODO: currently only applicable to 3-D cases
240 point p0 = lPosition;
241 scalar power(0.0);
242 scalar area(0.0);
243 p1 = p0;
244 if (mesh_.nGeometricD() == 3)
245 {
246 for (label ri = 0; ri < ndr_; ri++)
247 {
248 scalar r1 = SMALL + dr*ri;
249
250 scalar r2 = r1 + dr;
251
252 scalar rP = ((r1 + r2)/2);
253
254 // local radius on disk
255 vector localR = ((r1 + r2)/2)*rArea;
256
257 // local final radius on disk
258 vector finalR = rP*rArea;
259
260 scalar theta0 = 0.0;//dTheta/2.0;
261 for (label thetai = 0; thetai < ndTheta_; thetai++)
262 {
263 scalar theta1 = theta0 + SMALL + dTheta*thetai;
264
265 scalar theta2 = theta1 + dTheta;
266
267 scalar thetaP = (theta1 + theta2)/2.0;
268
269 quaternion Q(lDir, thetaP);
270
271 // Initial position on disk
272 vector initialPos = (Q.R() & localR);
273
274 // Final position on disk
275 vector finalPos = (Q.R() & finalR);
276
277 // Initial position
278 p0 = lPosition + initialPos;
279
280 // calculate target point using new deviation rl
281 p1 = lPosition + finalPos + (0.5*maxTrackLength_*lDir);
282
283 scalar Ip = calculateIp(rP, thetaP);
284
285 scalar dAi = (sqr(r2) - sqr(r1))*(theta2 - theta1)/2.0;
286
287 power += Ip*dAi;
288 area += dAi;
289
290 label cellI = mesh_.findCell(p0);
291
292 if (cellI != -1)
293 {
294 // Create a new particle
295 DTRMParticle* pPtr =
296 new DTRMParticle(mesh_, p0, p1, Ip, cellI, dAi, -1);
297
298 // Add to cloud
299 DTRMCloud_.addParticle(pPtr);
300 }
301
302 if (returnReduce(cellI, maxOp<label>()) == -1)
303 {
304 if (++nMissed <= 10)
305 {
307 << "Cannot find owner cell for focalPoint at "
308 << p0 << endl;
309 }
310 }
311 }
312 }
313 }
314 else
315 {
317 << "Current functionality limited to 3-D cases"
318 << exit(FatalError);
319 }
320
321 if (nMissed)
322 {
323 Info<< "Seeding missed " << nMissed << " locations" << endl;
324 }
325
327 << "Total Power in the laser : " << power << nl
328 << "Total Area in the laser : " << area << nl
329 << endl;
330}
331
332
333// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
334
336:
337 radiationModel(typeName, T),
338 mode_(powerDistNames_.get("mode", *this)),
339 DTRMCloud_(mesh_, "DTRMCloud", IDLList<DTRMParticle>()),
340 nParticles_(0),
341 ndTheta_(get<label>("nTheta")),
342 ndr_(get<label>("nr")),
343 maxTrackLength_(mesh_.bounds().mag()),
344
345 focalLaserPosition_
346 (
347 Function1<point>::New("focalLaserPosition", *this, &mesh_)
348 ),
349
350 laserDirection_
351 (
352 Function1<vector>::New("laserDirection", *this, &mesh_)
353 ),
354
355 focalLaserRadius_(get<scalar>("focalLaserRadius")),
356 qualityBeamLaser_
357 (
358 getOrDefault<scalar>("qualityBeamLaser", 0)
359 ),
360
361 sigma_(0),
362 I0_(0),
363 laserPower_(Function1<scalar>::New("laserPower", *this, &mesh_)),
364 powerDistribution_(),
365
366 reflectionSwitch_(false),
367
368 alphaCut_(getOrDefault<scalar>("alphaCut", 0.5)),
369
370 a_
371 (
373 (
374 "a",
375 mesh_.time().timeName(),
376 mesh_,
377 IOobject::NO_READ,
378 IOobject::NO_WRITE
379 ),
380 mesh_,
382 ),
383 e_
384 (
386 (
387 "e",
388 mesh_.time().timeName(),
389 mesh_,
390 IOobject::NO_READ,
391 IOobject::NO_WRITE
392 ),
393 mesh_,
395 ),
396 E_
397 (
399 (
400 "E",
401 mesh_.time().timeName(),
402 mesh_,
403 IOobject::NO_READ,
404 IOobject::NO_WRITE
405 ),
406 mesh_,
408 ),
409 Q_
410 (
412 (
413 "Q",
414 mesh_.time().timeName(),
415 mesh_,
416 IOobject::NO_READ,
417 IOobject::AUTO_WRITE
418 ),
419 mesh_,
421 )
422{
423 initialiseReflection();
424
425 initialise();
426}
427
428
430(
431 const dictionary& dict,
432 const volScalarField& T
433)
434:
435 radiationModel(typeName, dict, T),
436 mode_(powerDistNames_.get("mode", *this)),
437 DTRMCloud_(mesh_, "DTRMCloud", IDLList<DTRMParticle>()),
438 nParticles_(0),
439 ndTheta_(get<label>("nTheta")),
440 ndr_(get<label>("nr")),
441 maxTrackLength_(mesh_.bounds().mag()),
442
443 focalLaserPosition_
444 (
445 Function1<point>::New("focalLaserPosition", *this, &mesh_)
446 ),
447 laserDirection_
448 (
449 Function1<vector>::New("laserDirection", *this, &mesh_)
450 ),
451
452 focalLaserRadius_(get<scalar>("focalLaserRadius")),
453 qualityBeamLaser_
454 (
455 getOrDefault<scalar>("qualityBeamLaser", 0)
456 ),
457
458 sigma_(0),
459 I0_(0),
460 laserPower_(Function1<scalar>::New("laserPower", *this, &mesh_)),
461 powerDistribution_(),
462
463 reflectionSwitch_(false),
464
465 alphaCut_(getOrDefault<scalar>("alphaCut", 0.5)),
466
467 a_
468 (
470 (
471 "a",
472 mesh_.time().timeName(),
473 mesh_,
474 IOobject::NO_READ,
475 IOobject::NO_WRITE
476 ),
477 mesh_,
479 ),
480 e_
481 (
483 (
484 "e",
485 mesh_.time().timeName(),
486 mesh_,
487 IOobject::NO_READ,
488 IOobject::NO_WRITE
489 ),
490 mesh_,
492 ),
493 E_
494 (
496 (
497 "E",
498 mesh_.time().timeName(),
499 mesh_,
500 IOobject::NO_READ,
501 IOobject::NO_WRITE
502 ),
503 mesh_,
505 ),
506 Q_
507 (
509 (
510 "Q",
511 mesh_.time().timeName(),
512 mesh_,
513 IOobject::NO_READ,
514 IOobject::AUTO_WRITE
515 ),
516 mesh_,
518 )
519{
520 initialiseReflection();
521 initialise();
522}
523
524
525// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
526
528{
530 {
531 return true;
532 }
533
534 return false;
535}
536
537Foam::label Foam::radiation::laserDTRM::nBands() const
538{
539 return 1;
540}
541
542
544{
545 tmp<volScalarField> treflectingCells
546 (
548 (
550 (
551 "reflectingCellsVol",
552 mesh_.time().timeName(),
553 mesh_,
556 ),
557 mesh_,
558 dimensionedScalar("zero", dimless, -1)
559 )
560 );
561 volScalarField& reflectingCellsVol = treflectingCells.ref();
562
563
565 (
567 (
569 (
570 "nHat",
571 mesh_.time().timeName(),
572 mesh_,
575 ),
576 mesh_,
578 )
579 );
580 volVectorField& nHat = tnHat.ref();
581
582
583 // Reset the field
584 Q_ == dimensionedScalar(Q_.dimensions(), Zero);
585
586 a_ = absorptionEmission_->a();
587 e_ = absorptionEmission_->e();
588 E_ = absorptionEmission_->E();
589
590 const interpolationCell<scalar> aInterp(a_);
591 const interpolationCell<scalar> eInterp(e_);
592 const interpolationCell<scalar> EInterp(E_);
593 const interpolationCell<scalar> TInterp(T_);
594
595 labelField reflectingCells(mesh_.nCells(), -1);
596
597 UPtrList<reflectionModel> reflectionUPtr;
598
599 if (reflectionSwitch_)
600 {
601 reflectionUPtr.resize(reflections_.size());
602
603 label reflectionModelId(0);
604 forAllIters(reflections_, iter1)
605 {
606 reflectionModel& model = iter1()();
607
608 reflectionUPtr.set(reflectionModelId, &model);
609
610 const volScalarField& alphaFrom =
611 mesh_.lookupObject<volScalarField>
612 (
613 IOobject::groupName("alpha", iter1.key().first())
614 );
615
616 const volScalarField& alphaTo =
617 mesh_.lookupObject<volScalarField>
618 (
619 IOobject::groupName("alpha", iter1.key().second())
620 );
621
622 const volVectorField nHatPhase(nHatfv(alphaFrom, alphaTo));
623
624 const volScalarField gradAlphaf
625 (
626 fvc::grad(alphaFrom)
627 & fvc::grad(alphaTo)
628 );
629
630 const volScalarField nearInterface(pos(alphaTo - alphaCut_));
631
632 const volScalarField mask(nearInterface*gradAlphaf);
633
634 forAll(alphaFrom, cellI)
635 {
636 if
637 (
638 nearInterface[cellI]
639 && mag(nHatPhase[cellI]) > 0.99
640 && mask[cellI] < 0
641 )
642 {
643 reflectingCells[cellI] = reflectionModelId;
644 reflectingCellsVol[cellI] = reflectionModelId;
645 if (mag(nHat[cellI]) == 0.0)
646 {
647 nHat[cellI] += nHatPhase[cellI];
648 }
649 }
650 }
651 reflectionModelId++;
652 }
653 }
654
655 interpolationCellPoint<vector> nHatInterp(nHat);
656
658 (
659 DTRMCloud_,
660 aInterp,
661 eInterp,
662 EInterp,
663 TInterp,
664 nHatInterp,
665 reflectingCells,
666 reflectionUPtr,
667 Q_
668 );
669
670 Info<< "Move particles..."
671 << returnReduce(DTRMCloud_.size(), sumOp<label>()) << endl;
672
673 DTRMCloud_.move(DTRMCloud_, td, mesh_.time().deltaTValue());
674
675 // Normalize by cell volume
676 Q_.primitiveFieldRef() /= mesh_.V();
677
678 if (debug)
679 {
680 Info<< "Final number of particles..."
681 << returnReduce(DTRMCloud_.size(), sumOp<label>()) << endl;
682
683 pointField lines(2*DTRMCloud_.size());
684 {
685 label i = 0;
686 for (const DTRMParticle& p : DTRMCloud_)
687 {
688 lines[i] = p.position();
689 lines[i+1] = p.p0();
690 i += 2;
691 }
692 }
693
695
696 if (Pstream::master())
697 {
698 OBJstream os(type() + ":particlePath.obj");
699
700 for (label pointi = 0; pointi < lines.size(); pointi += 2)
701 {
702 os.write(linePointRef(lines[pointi], lines[pointi+1]));
703 }
704 }
705
706 scalar totalQ = gSum(Q_.primitiveFieldRef()*mesh_.V());
707 Info << "Total energy absorbed [W]: " << totalQ << endl;
708
709 if (mesh_.time().outputTime())
710 {
711 reflectingCellsVol.write();
712 nHat.write();
713 }
714 }
715
716 // Clear and initialise the cloud
717 // NOTE: Possible to reset original particles, but this requires
718 // data transfer for the cloud in differet processors.
719 initialise();
720}
721
722
724{
726 (
728 (
729 "zero",
730 mesh_.time().timeName(),
731 mesh_,
734 false
735 ),
736 mesh_,
738 );
739}
740
741
744{
745 return Q_.internalField();
746}
747
748
749// ************************************************************************* //
bool found
Macros for easy insertion into run-time selection tables.
const volScalarField & alpha1
const volScalarField & alpha2
Base cloud calls templated on particle type.
Definition: Cloud.H:68
Class used to pass tracking data to the trackToFace function.
Definition: DTRMParticle.H:95
Discrete Transfer Radiation Model (DTRM) particle.
Definition: DTRMParticle.H:65
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: Function1.H:96
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
Template class for intrusive linked lists.
Definition: ILList.H:69
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
OFstream that keeps track of vertices.
Definition: OBJstream.H:61
Random number generator.
Definition: Random.H:60
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
Definition: UPtrList.H:71
const T * set(const label i) const
Definition: UPtrList.H:248
void resize(const label newLen)
Change the size of the list.
Definition: UPtrListI.H:190
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
const Type & value() const
Return const reference to value.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
static void gatherInplaceOp(List< Type > &fld, const int tag=UPstream::msgType(), const UPstream::commsTypes=UPstream::commsTypes::nonBlocking, const label comm=UPstream::worldComm)
Inplace collect data in processor order on master (in serial: a no-op).
2D table interpolation. The data must be in ascending order in both dimensions x and y.
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
Uses the cell value for any location within the cell.
An ordered or unorder pair of phase names. Typically specified as follows.
Definition: phasePairKey.H:68
Quaternion class used to perform rotations in 3D space.
Definition: quaternion.H:58
Discrete Tray Radiation Method for collimated radiation flux. At the moment the particles are injecte...
Definition: laserDTRM.H:81
virtual label nBands() const
Number of bands for this radiation model.
virtual tmp< DimensionedField< scalar, volMesh > > Ru() const
Source term component (constant)
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
bool read()
Read radiation properties dictionary.
void calculate()
Solve radiation equation(s)
Lookup type of boundary radiation properties.
Definition: lookup.H:66
Top level model for radiation modelling.
const fvMesh & mesh_
Reference to the mesh database.
virtual bool read()=0
Read radiationProperties dictionary.
Base class for radiation scattering.
virtual bool write(const bool valid=true) const
Write using setting from DB.
splitCell * master() const
Definition: splitCell.H:113
A class for managing temporary objects.
Definition: tmp.H:65
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define defineTemplateTypeNameAndDebugWithName(Type, Name, DebugSwitch)
Define the typeName and debug information, lookup as Name.
Definition: className.H:126
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
volScalarField & p
const volScalarField & T
const volScalarField & p0
Definition: EEqn.H:36
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
Calculate the matrix for the laplacian of the field.
Calculate the finiteVolume matrix for implicit and explicit sources.
word timeName
Definition: getTimeIndex.H:3
#define DebugInfo
Report an information message using Foam::Info.
#define WarningInFunction
Report a warning using Foam::Warning.
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
Different types of constants.
const wordList area
Standard area field types (scalar, vector, tensor, etc)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
auto key(const Type &t) -> typename std::enable_if< std::is_enum< Type >::value, typename std::underlying_type< Type >::type >::type
Definition: foamGltfBase.H:108
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet dimPower
Type gSum(const FieldField< Field, Type > &f)
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:680
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
line< point, const point & > linePointRef
A line using referred points.
Definition: linePointRef.H:47
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:54
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensionedScalar pow4(const dimensionedScalar &ds)
errorManip< error > abort(error &err)
Definition: errorManip.H:144
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
error FatalError
dimensionedScalar cbrt(const dimensionedScalar &ds)
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
#define addToRadiationRunTimeSelectionTables(model)
dictionary dict
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllIters(container, iter)
Iterate across all elements in the container object.
Definition: stdFoam.H:260
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278
A non-counting (dummy) refCount.
Definition: refCount.H:59
Unit conversion functions.