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 // Populate solid field pointer lists
11 {
12 Info<< "*** Reading solid mesh thermophysical properties for region "
13 << solidRegions[i].name() << nl << endl;
14
15 Info<< " Adding to thermos\n" << endl;
16 thermos.set(i, solidThermo::New(solidRegions[i]));
17
18 Info<< " Adding to radiations\n" << endl;
19 radiations.set(i, radiation::radiationModel::New(thermos[i].T()));
20
21 Info<< " Adding fvOptions\n" << endl;
23 (
24 i,
25 new fv::options(solidRegions[i])
26 );
27
28 if (!thermos[i].isotropic())
29 {
30 Info<< " Adding coordinateSystems\n" << endl;
32 (
33 i,
34 coordinateSystem::New
35 (
36 solidRegions[i],
37 thermos[i],
38 coordinateSystem::typeName_()
39 )
40 );
41
42 tmp<volVectorField> tkappaByCp =
43 thermos[i].Kappa()/thermos[i].Cp();
44
45 aniAlphas.set
46 (
47 i,
48 new volSymmTensorField
49 (
50 IOobject
51 (
52 "Anialpha",
53 runTime.timeName(),
54 solidRegions[i],
55 IOobject::NO_READ,
56 IOobject::NO_WRITE
57 ),
58 solidRegions[i],
59 dimensionedSymmTensor(tkappaByCp().dimensions(), Zero),
60 zeroGradientFvPatchSymmTensorField::typeName
61 )
62 );
63
64 aniAlphas[i].primitiveFieldRef() =
65 coordinates[i].transformPrincipal
66 (
67 solidRegions[i].cellCentres(),
68 tkappaByCp()
69 );
70 aniAlphas[i].correctBoundaryConditions();
71
72 }
73
74 IOobject betavSolidIO
75 (
76 "betavSolid",
77 runTime.timeName(),
78 solidRegions[i],
79 IOobject::MUST_READ,
80 IOobject::AUTO_WRITE
81 );
82
83 if (betavSolidIO.typeHeaderOk<volScalarField>(true))
84 {
85 betavSolid.set
86 (
87 i,
88 new volScalarField(betavSolidIO, solidRegions[i])
89 );
90 }
91 else
92 {
93 betavSolid.set
94 (
95 i,
96 new volScalarField
97 (
98 IOobject
99 (
100 "betavSolid",
101 runTime.timeName(),
102 solidRegions[i],
103 IOobject::NO_READ,
104 IOobject::NO_WRITE
105 ),
106 solidRegions[i],
107 dimensionedScalar("1", dimless, scalar(1))
108 )
109 );
110 }
111 }
const T * set(const label i) const
Definition: PtrList.H:138
const volScalarField & T
engineTime & runTime
PtrList< coordinateSystem > coordinates(solidRegions.size())
PtrList< solidThermo > thermos(solidRegions.size())
PtrList< volSymmTensorField > aniAlphas(solidRegions.size())
PtrList< volScalarField > betavSolid(solidRegions.size())
PtrList< fv::options > solidHeatSources(solidRegions.size())
PtrList< radiation::radiationModel > radiations(solidRegions.size())
PtrList< fvMesh > solidRegions(solidNames.size())
IOobject betavSolidIO("betavSolid", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333