StCorr.H
Go to the documentation of this file.
1  dimensionedScalar StCorr("StCorr", dimless, 1.0);
2 
3  if (ign.igniting())
4  {
5  // Calculate volume of ignition kernel
6  dimensionedScalar Vk("Vk", dimVolume, gSum(c*mesh.V().field()));
7  dimensionedScalar Ak("Ak", dimArea, 0.0);
8 
9  if (Vk.value() > small)
10  {
11  // Calculate kernel area from its volume
12  // and the dimensionality of the case
13 
14  switch(mesh.nGeometricD())
15  {
16  case 3:
17  {
18  // Assume it is part-spherical
19  scalar sphereFraction
20  (
22  (
23  combustionProperties.lookup
24  (
25  "ignitionSphereFraction"
26  )
27  )
28  );
29 
30  Ak = sphereFraction*4.0*constant::mathematical::pi
31  *pow
32  (
33  3.0*Vk
34  /(sphereFraction*4.0*constant::mathematical::pi),
35  2.0/3.0
36  );
37  }
38  break;
39 
40  case 2:
41  {
42  // Assume it is part-circular
43  dimensionedScalar thickness
44  (
45  combustionProperties.lookup("ignitionThickness")
46  );
47 
48  scalar circleFraction
49  (
51  (
52  combustionProperties.lookup
53  (
54  "ignitionCircleFraction"
55  )
56  )
57  );
58 
59  Ak = circleFraction*constant::mathematical::pi*thickness
60  *sqrt
61  (
62  4.0*Vk
63  /(
64  circleFraction
65  *thickness
67  )
68  );
69  }
70  break;
71 
72  case 1:
73  // Assume it is plane or two planes
75  (
76  combustionProperties.lookup("ignitionKernelArea")
77  );
78  break;
79  }
80 
81  // Calculate kernel area from b field consistent with the
82  // discretisation of the b equation.
83  const volScalarField mgb
84  (
85  fvc::div(nf, b, "div(phiSt,b)") - b*fvc::div(nf) + dMgb
86  );
87  dimensionedScalar AkEst = gSum(mgb*mesh.V().field());
88 
89  StCorr.value() = max(min((Ak/AkEst).value(), 10.0), 1.0);
90 
91  Info<< "StCorr = " << StCorr.value() << endl;
92  }
93  }
#define readScalar
Definition: doubleScalar.C:38
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
dimensionedScalar StCorr("StCorr", dimless, 1.0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
dynamicFvMesh & mesh
Type gSum(const FieldField< Field, Type > &f)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionedScalar c
Speed of light in a vacuum.
messageStream Info
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57