# Thermal Stability Investigation in a Reactive Sphere of Combustible Material

Authors: R. S. Lebelo
Journal: Advances in Mathematical Physics
 (2016)

---

## Abstract

An investigation of thermal stability in a stockpile of combustible material is considered. The combustible material is any carbon containing material that can react with oxygen trapped in a stockpile due to exothermic chemical reaction. The complicated process is modelled in a sphere and one-dimensional energy equation is used to solve the problem. The semi-implicit finite difference method (FDM) is applied to tackle the nonlinear differential equation governing the problem. Graphical solutions are displayed to describe effects of embedded kinetic parameters on the temperature of the system.


---

## Body

## 1. Introduction

Thermal stability in a stockpile of combustible material has drawn attention of many researchers. The interest is sparked by self-ignited fires that cause hazards to the environment, vegetation, and domestic and wild animals. Some veld fires are caused by stockpiles of combustible materials that have been left without attention for a long time. Some physical factors or parameters that influence the self-ignition process are combustible material particle size, volume-to-surface ratio, porosity, thermal conductivity, density, and heat capacity including convection in the surroundings [1–3]. The three groups, into which the parameters that influence self-ignition in stockpiles of combustible materials are distinguished, are outlined in [4]. Ignition also depends on chemical kinetic factors, and in this paper we study the effects of these embedded parameters on thermal stability of self-ignited processes. Exothermic chemical reaction taking place in a reactive stockpile results in heat generation. Should the heat generation within a stockpile exceed heat loss to the ambient surroundings, thermal instability and runaway leading to self-ignition may occur [5–7]. Thermal explosion criticality concept was investigated by Frank-Kamenetskii and he also developed the steady state theory for combustion of reactive materials due to exothermic chemical reaction [8]. Frank-Kamenetskii’s work enables researchers to evaluate thermal criticality values for steady state exothermic chemical reaction. These criticality values are helpful to study thermal stability in processes of self-ignition of combustible materials due to exothermic chemical reaction. Self-ignited combustion is a complicated process that involves many radicals which are short-lived and their interaction is nonlinear [9, 10]. The complicated processes of combustion due to exothermic chemical reactions in reactive stockpile are modelled mathematically in order to assess materials property thermal stability, climate change indicators, and mitigation strategies [11, 12]. Thermal stability of reactive materials was studied by several authors [13–15] who carried out the investigation in a reactive slab. In this paper, the investigation is carried out in a sphere of reactive material. The nonlinear differential equation governing the problem is tackled numerically using semi-implicit FDM. The paper layout is done as follows; in Section 2, the mathematical formulation is described. The numerical analysis is done in Section 3 and results and discussion are given in Section 4.


## 2. Mathematical Formulation

The study of thermal stability of a stockpile of reactive material is modelled in a sphere of constant thermal conductivity. It is assumed that the sphere undergoes annth-order oxidation chemical reaction and that no reaction consumption takes place. Figure 1 illustrates the geometry of the problem.Figure 1
Geometry of the problem.The complicated chemistry involved in this process is simplified by assuming one-step finite rate irreversible Arrhenius kinetics. The one-dimensional nonlinear partial differential equation that governs the problem is [3](1)
p
c
p
∂
T
∂
t
-
=
k
1
r
-
2
∂
∂
r
-
r
-
2
∂
T
∂
r
+
Q
A
C
K
T
v
l
m
e
-
E
/
R
T
-
φ
T
-
T
b
.The initial condition is(2)
T
r
-
,
0
=
T
0
,with boundary conditions(3)
∂
T
∂
r
-
0
,
t
-
=
0
;
T
a
,
t
-
=
T
b
.In this case, T is the sphere’s absolute temperature, T
0 is the initial temperature of the sphere, and T
b is the ambient temperature. p is the density, c
p is the specific heat at constant pressure, k is the thermal conductivity of the material, Q is the heat of reaction, A is the rate constant, and C is the reactant concentration. We have also K as Boltzmann’s constant, v as the vibration frequency, l as Planck’s number, E as the activation energy, and R as the universal gas number; m is the numerical exponent and takes the following values: −2 for sensitized (light induced), 0 for Arrhenius, and 1
/
2 for bimolecular kinetics; φ is the heat loss parameter [3, 8, 9].The following dimensionless parameters are introduced to (1)–(3): (4)
θ
=
E
T
-
T
b
R
T
b
2
,
θ
0
=
E
T
0
-
T
b
R
T
b
2
,
r
=
r
-
a
,
ε
=
R
T
b
E
,
t
=
k
t
-
p
c
p
a
2
,
δ
=
a
2
φ
k
,
λ
=
K
T
b
v
l
m
Q
A
E
a
2
C
k
R
T
b
2
exp
⁡
-
E
R
T
b
.Equations (1)–(3) take the following form:(5)
∂
θ
∂
t
=
1
r
2
∂
∂
r
r
2
∂
θ
∂
r
+
λ
1
+
ε
θ
m
e
θ
/
1
+
ε
θ
-
δ
θ
.The initial condition is(6)
θ
r
,
0
=
θ
0
,and boundary conditions are(7)
∂
θ
∂
r
0
,
t
=
0
;
θ
1
,
t
=
0
,where θ is the dimensionless temperature, θ
0 is the dimensionless initial temperature, λ is the Frank-Kamenetskii parameter (reaction rate parameter), ε is the activation energy parameter, r is the dimensionless radial distance, and δ is the heat loss parameter.


## 3. Numerical Approach

Following [5], the semi-implicit finite difference scheme is applied where implicit terms are taken at intermediate time level (
N
+
ξ
) for 0
≤
ξ
≤
1. In this investigation, we use ξ
=
1 to enable use of larger time steps that allow working with any value of the time step. The governing partial differential equation is approximated with both the second and first spatial derivatives with second-order central differences. Boundary conditions are obtained from modification of equations corresponding to the first and last grid points. The governing equation’s discretization is based on linear Cartesian mesh and uniform grid on which finite differences are taken.Each component of the governing PDE, wherej represents position, N time, and Δ
r mesh spacing, is transformed as follows:(8)
∂
θ
∂
t
=
θ
j
N
+
1
-
θ
j
N
Δ
t
,
∂
θ
∂
r
=
1
2
Δ
r
θ
j
+
1
N
+
1
-
θ
j
-
1
N
+
1
+
θ
j
+
1
N
-
θ
j
-
1
N
,
∂
2
θ
∂
r
2
=
1
Δ
r
2
ξ
θ
j
+
1
N
+
1
-
2
θ
j
N
+
1
+
θ
j
-
1
N
+
1
+
1
-
ξ
θ
j
+
1
N
-
2
θ
j
N
+
θ
j
-
1
N
.It follows that the semi-implicit FDM for the governing equation is expressed as follows:(9)
θ
N
+
1
-
θ
N
Δ
t
=
∂
2
∂
r
2
θ
N
+
ξ
+
2
r
∂
∂
r
θ
N
+
λ
1
+
ε
θ
m
e
θ
/
1
+
ε
θ
N
-
δ
θ
N
.The new time terms N
+
1 are put on the left. These are the unknowns. The known time terms, that is, the present time terms N, are put on the right to get the following:(10)
-
ξ
τ
θ
j
+
1
N
+
1
+
1
+
2
ξ
τ
θ
j
N
+
1
-
ξ
τ
θ
j
-
1
N
+
1
=
τ
1
-
ξ
θ
j
+
1
N
+
1
-
2
τ
1
-
ξ
-
δ
θ
j
N
+
τ
1
-
ξ
θ
j
-
1
N
+
1
2
r
j
2
τ
θ
j
+
1
N
-
θ
j
-
1
N
+
λ
1
+
ε
θ
m
e
θ
/
1
+
ε
θ
N
,where τ
=
Δ
t
/
Δ
r
2. Equation (10) forms a tridiagonal system. Maple software was applied in this study to give graphical solutions.


## 4. Results and Discussion

This section clarifies understanding of thermal development in the combustion process due to exothermic chemical reaction. Numerical values were assigned to embedded kinetic parameters to investigate their effect on the temperature of the system. It is assumed that, att
=
0, the sphere’s temperature is equal to that of the ambient surroundings and that there is no convective heat loss to the surroundings. But the temperature of the system increases at t
>
0 due to progressive exothermic chemical reaction within the system that results with heat loss from the sphere since its temperature will be higher than that of the surroundings.


### 4.1. Effects of Kinetic Parameters on Temperature

We see from Figures2–5 the behavior of temperature in response to variation of kinetic parameters embedded in the system. Figure 2 illustrates the behavior of temperature as t increases while parameters, λ
,
ε
,
m, and δ, are kept constant. Generally, as heat loss from the sphere to the ambient surroundings takes place, the temperature is highest at the center and lowest on the surface. From Figure 2, we observe also that the sphere’s temperature attains its steady state value. Once the steady state value is reached, there is no more increase in temperature of the system. We observe in Figure 3 how temperature varies with λ, the rate of reaction parameter. An increase in λ results with corresponding increase of temperature fields of the sphere due to progressive exothermic chemical reaction of the system. Therefore, high values of temperature are attained, and, should λ keep on increasing without corresponding heat loss to the surroundings, the increasing accumulated heat within the system may lead to thermal runaways to a point where ignition occurs. The same scenario is observed in Figure 4 as ε increases. We observe also that ε has little effect on temperature behavior compared to λ. Use of smaller values to increase ε does not clearly show its effect on temperature behavior, hence the use of higher values. The reason of this behavior is that the exothermic chemical reaction is assumed to take place in a spherical domain which, in general, is known to have the highest thermal stability as compared to cylindrical and rectangular domains. In other words, the transfer of heat from hot surroundings to a liquid inside a spherical vessel is lesser than that for a cylinder or a rectangular container. This is so because a spherical vessel has smaller surface area per unit volume than cylindrical and rectangular vessels.Figure 2
Temperature with time.Figure 3
Temperature withλ.Figure 4
Temperature withε.Figure 5
Temperature withδ.An increase inδ, the heat loss parameter, gives a different scenario as illustrated by Figure 5, which shows a decrease in temperature fields. This parameter affirms that considerable loss of heat to the surroundings enhances thermal stability and therefore minimizes thermal runaways to reduce ignition.


### 4.2. Thermal Stability and Blow-Ups

In this section, we look at plots of Nusselt number (N
u) versus Frank Kamenetskii parameter λ, for varied values of δ
,
ε, and m, respectively. We see from Figures 6 and 7 that thermal stability is attained by keeping high values of δ and ε. A different scenario is illustrated by Figure 8. In this case, thermal stability is attained by keeping m very low. Table 1 provides numerical values obtained for N
u and λ, for different values of δ
,
ε, and m, respectively.Table 1
Numerical values showing the effects ofδ, ε, and m on thermal criticality values.
δ
ε
m
N
u
λ
c
1
0.1
0.5
2.39524
3.9197
3
0.1
0.5
2.46169
4.6458
5
0.1
0.5
2.53119
5.3734
1
0.5
0.5
14.70735
7.6589
1
0.8
0.5
7.75651
7.6493
1
0.1
0
2.56740
4.1587
1
0.1
−2
3.70737
5.5635Figure 6
δ on thermal stability.Figure 7
ε on thermal stability.Figure 8
m on thermal stability.


## 4.1. Effects of Kinetic Parameters on Temperature

We see from Figures2–5 the behavior of temperature in response to variation of kinetic parameters embedded in the system. Figure 2 illustrates the behavior of temperature as t increases while parameters, λ
,
ε
,
m, and δ, are kept constant. Generally, as heat loss from the sphere to the ambient surroundings takes place, the temperature is highest at the center and lowest on the surface. From Figure 2, we observe also that the sphere’s temperature attains its steady state value. Once the steady state value is reached, there is no more increase in temperature of the system. We observe in Figure 3 how temperature varies with λ, the rate of reaction parameter. An increase in λ results with corresponding increase of temperature fields of the sphere due to progressive exothermic chemical reaction of the system. Therefore, high values of temperature are attained, and, should λ keep on increasing without corresponding heat loss to the surroundings, the increasing accumulated heat within the system may lead to thermal runaways to a point where ignition occurs. The same scenario is observed in Figure 4 as ε increases. We observe also that ε has little effect on temperature behavior compared to λ. Use of smaller values to increase ε does not clearly show its effect on temperature behavior, hence the use of higher values. The reason of this behavior is that the exothermic chemical reaction is assumed to take place in a spherical domain which, in general, is known to have the highest thermal stability as compared to cylindrical and rectangular domains. In other words, the transfer of heat from hot surroundings to a liquid inside a spherical vessel is lesser than that for a cylinder or a rectangular container. This is so because a spherical vessel has smaller surface area per unit volume than cylindrical and rectangular vessels.Figure 2
Temperature with time.Figure 3
Temperature withλ.Figure 4
Temperature withε.Figure 5
Temperature withδ.An increase inδ, the heat loss parameter, gives a different scenario as illustrated by Figure 5, which shows a decrease in temperature fields. This parameter affirms that considerable loss of heat to the surroundings enhances thermal stability and therefore minimizes thermal runaways to reduce ignition.


## 4.2. Thermal Stability and Blow-Ups

In this section, we look at plots of Nusselt number (N
u) versus Frank Kamenetskii parameter λ, for varied values of δ
,
ε, and m, respectively. We see from Figures 6 and 7 that thermal stability is attained by keeping high values of δ and ε. A different scenario is illustrated by Figure 8. In this case, thermal stability is attained by keeping m very low. Table 1 provides numerical values obtained for N
u and λ, for different values of δ
,
ε, and m, respectively.Table 1
Numerical values showing the effects ofδ, ε, and m on thermal criticality values.
δ
ε
m
N
u
λ
c
1
0.1
0.5
2.39524
3.9197
3
0.1
0.5
2.46169
4.6458
5
0.1
0.5
2.53119
5.3734
1
0.5
0.5
14.70735
7.6589
1
0.8
0.5
7.75651
7.6493
1
0.1
0
2.56740
4.1587
1
0.1
−2
3.70737
5.5635Figure 6
δ on thermal stability.Figure 7
ε on thermal stability.Figure 8
m on thermal stability.


## 5. Conclusion

In this paper, we investigated thermal stability of a reactive sphere of combustible material. It is assumed that the sphere is of constant thermal conductivity and that no reactant consumption takes place. The kinetic parametersλ and ε were found to enhance exothermic chemical reaction of the system and therefore reduce thermal stability. Contrary to these parameters, δ retards the exothermic chemical reaction and hence boosts thermal stability. Thermal stability is also achieved by using the lowest values of m and high values of δ and ε. Theoretical investigation of thermal stability in stockpiles of combustible materials provides a cost-effective, easier, and faster means than the experimental approach. The limitation to the model in this investigation is that only theoretical approach was considered, and therefore physical parameters mentioned in Section 1 cannot be dealt with. This investigation can be extended to systems where reactant consumption takes place and the sphere is of variable thermal conductivity.



---

*Source: 9384541-2016-07-31.xml*