createSolidFields.H
Go to the documentation of this file.
1  // Initialise solid field pointer lists
2  PtrList<coordinateSystem> coordinates(solidRegions.size());
3  PtrList<solidThermo> thermos(solidRegions.size());
4  PtrList<radiation::radiationModel> radiations(solidRegions.size());
5  PtrList<fv::options> solidHeatSources(solidRegions.size());
6  PtrList<volScalarField> betavSolid(solidRegions.size());
7  PtrList<volSymmTensorField> aniAlphas(solidRegions.size());
8 
9  List<bool> residualReachedSolid(solidRegions.size(), true);
10  List<bool> residualControlUsedSolid(solidRegions.size(), false);
11 
12  // Populate solid field pointer lists
14  {
15  Info<< "*** Reading solid mesh thermophysical properties for region "
16  << solidRegions[i].name() << nl << endl;
17 
18  Info<< " Adding to thermos\n" << endl;
20 
21  Info<< " Adding to radiations\n" << endl;
23 
24  Info<< " Adding fvOptions\n" << endl;
26  (
27  i,
28  new fv::options(solidRegions[i])
29  );
30 
31  if (!thermos[i].isotropic())
32  {
33  Info<< " Adding coordinateSystems\n" << endl;
35  (
36  i,
38  );
39 
40  tmp<volVectorField> tkappaByCp =
41  thermos[i].Kappa()/thermos[i].Cp();
42 
43  aniAlphas.set
44  (
45  i,
47  (
48  IOobject
49  (
50  "Anialpha",
51  runTime.timeName(),
52  solidRegions[i],
53  IOobject::NO_READ,
54  IOobject::NO_WRITE
55  ),
56  solidRegions[i],
58  (
59  "zero",
60  tkappaByCp().dimensions(),
61  Zero
62  ),
63  zeroGradientFvPatchSymmTensorField::typeName
64  )
65  );
66 
67  aniAlphas[i].primitiveFieldRef() =
68  coordinates[i].R().transformVector(tkappaByCp());
69  aniAlphas[i].correctBoundaryConditions();
70 
71  }
72 
73  IOobject betavSolidIO
74  (
75  "betavSolid",
76  runTime.timeName(),
77  solidRegions[i],
78  IOobject::MUST_READ,
79  IOobject::AUTO_WRITE
80  );
81 
82  if (betavSolidIO.typeHeaderOk<volScalarField>(true))
83  {
84  betavSolid.set
85  (
86  i,
87  new volScalarField(betavSolidIO, solidRegions[i])
88  );
89  }
90  else
91  {
92  betavSolid.set
93  (
94  i,
95  new volScalarField
96  (
97  IOobject
98  (
99  "betavSolid",
100  runTime.timeName(),
101  solidRegions[i],
102  IOobject::NO_READ,
103  IOobject::NO_WRITE
104  ),
105  solidRegions[i],
106  dimensionedScalar("1", dimless, scalar(1))
107  )
108  );
109  }
110  }
PtrList< volScalarField > betavSolid(solidRegions.size())
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:58
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
List< bool > residualReachedSolid(solidRegions.size(), true)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
List< bool > residualControlUsedSolid(solidRegions.size(), false)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
PtrList< fv::options > solidHeatSources(solidRegions.size())
PtrList< solidThermo > thermos(solidRegions.size())
static const zero Zero
Definition: zero.H:91
PtrList< volSymmTensorField > aniAlphas(solidRegions.size())
static const char nl
Definition: Ostream.H:262
dimensioned< symmTensor > dimensionedSymmTensor
Dimensioned tensor obtained from generic dimensioned type.
const volScalarField & T
forAll(solidRegions, i)
PtrList< coordinateSystem > coordinates(solidRegions.size())
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< fvMesh > solidRegions(solidsNames.size())
messageStream Info
PtrList< radiation::radiationModel > radiations(solidRegions.size())