# Kinetic development of low-temperature propane oxidation in a repetitively-
pulsed nanosecond discharge

Zhenyang Li  Bo Yin  Qifu Lin  Yifei Zhu<EMAIL_ADDRESS>Yun Wu

###### Abstract

The kinetics of plasma assisted low temperature oxidation of $\rm
C_{3}H_{8}/O_{2}/Ar$ mixtures have been studied in a wide specific deposition
energy with the help of nanosecond repetitively pulsed discharge. Two types of
nanosecond pulsed plasma sources, the nanosecond capillary discharge (nCD) and
dielectric barrier discharge (DBD) combined with the synchrotron
photoionization mass spectrometry are investigated. The electron impact
reaction rate of propane dissociation and some combustion chemical reactions
rate constants are updated according to the nCD and DBD experiment results,
and uncertainty of the reactions are analyzed in detail. Compared to the
existing model, the updated model’s prediction accuracy has great improvement
in species $\rm H_{2}O$, $\rm CO$, $\rm CO_{2}$, $\rm CH_{4}$, $\rm CH_{2}O$,
$\rm CH_{3}OH$, $\rm C_{2}H_{2}$, $\rm C_{2}H_{4}$, $\rm C_{2}H_{6}$, $\rm
C_{2}H_{5}OH$, $\rm C_{2}H_{5}OOH$, $\rm C_{3}H_{4}$-$\rm A$, $\rm
C_{3}H_{4}$-$\rm P$, $\rm C_{2}H_{5}CHO$, $\rm i$-$\rm C_{3}H_{7}OH$ and $\rm
C_{3}H_{7}OOH$. The propane oxidation processes assisted by DBD and nCD were
compared under different single pulse deposition energy (SPDE) conditions
while maintaining the same total deposition energy. The reduced electric field
in nCD is concentrated at 150-200 Td and 450-500 Td, whereas in DBD it ranges
from 0-25 Td and 50-250 Td. Notably, SPDE shows minimal influence on the
propane oxidation process, which is primarily controlled by total deposition
energy and demonstrates little dependence on the discharge type (DBD or nCD).

###### keywords:

Plasma-assisted combustion , Low-temperature oxidation , Propane , Kinetic
modeling

[1]organization=Institute of Energy, Hefei Comprehensive National Science
Center, city=Hefei, postcode=230026, country=China [2]organization=National
Key Lab of Aerospace Power System and Plasma Technology, Xi’an Jiaotong
University, city=Xi’an, postcode=710049, country=China
[3]organization=National Key Lab of Aerospace Power System and Plasma
Technology, Airforce Engineering University, city=Xi’an, postcode=710038,
country=China

## Novelty and significance statement

The kinetic of plasma assisted $\rm C_{3}H_{8}/O_{2}/Ar$ oxidation is
developed and studied in a wide specific deposition energy of nanosecond
repetitively pulsed discharge. Two types of nanosecond pulsed plasma sources
(nCD and DBD) are used to validate the kinetic mechanism. The electron impact
reaction rate of propane dissociation and some combustion chemical reactions
rate constants are updated according to the nCD and DBD experiment results.
The prediction accuracy has great improvement, and uncertainty of the
reactions are analyzed in detail. The propane oxidation processes assisted by
DBD and nCD were compared, and the propane oxidation processes is primarily
controlled by total deposition energy and demonstrates little dependence on
the discharge type (DBD or nCD).

## Authors contributions

Zhenyang Li: Writing-original draft, Formal analysis, Data curation. Bo Yin:
Investigation, Experiment development. Qifu Lin: Investigation, Supervision.
Yifei Zhu: Writing-review & editing, Conceptualization. Yun Wu: Funding
acquisition, Project administration.

##  1 Introduction

Non-equilibrium plasma-assisted combustion (PAC) offers an effective approach
for controlling ignition in internal combustion engines, industrial burners,
and aviation engines, thereby enhancing flame stability, reducing pollutant
emissions, and expanding ignition limits [1, 2, 3]. Understanding the detailed
mechanism of PAC requires studying the kinetic pathways of PAC and quantifying
their reaction rates. The kinetics of non-equilibrium plasma-assisted small
molecules like methane have been extensively studied [4, 5, 6, 7, 8].

Propane, as the smallest n-alkane with low-temperature reactivity, serves as a
typical molecule for studying the kinetics of large alkanes [9, 10],
attracting significant attention from researchers. Adamovich et al. [11]
proposed a plasma-assisted propane-air combustion kinetic model, which showed
good agreement with experimental temperature and time-resolved OH density
measurements. They further studied the plasma-assisted propane oxidation
kinetics in $\rm C_{3}H_{8}/O_{2}/Ar$ mixtures [12], validated the kinetic
model with experiments on alkane small molecule density at different
temperatures, and found large discrepancies for $\rm C_{2}H_{2}$ and $\rm
CH_{3}CHO$. Chen et al. [13] employed molecular beam mass spectrometry with
tunable synchrotron vacuum ultraviolet photoionization (SVUV-PIMS) to measure
hydrocarbons and oxygenated intermediates during plasma-assisted $\rm
C_{3}H_{8}/O_{2}/Ar$ processes, reporting numerous oxygenated intermediates
and establishing a kinetic model. However, significant discrepancies existed
between the model predictions and experimental measurements for certain
hydrocarbons and intermediates. Ban et al. [14] developed a propane-air plasma
combustion mechanism and validated it against OH density measurements from
experiments [15] that lacked sufficient validation. Based on the above
analysis, current kinetic models still exhibit significant discrepancies with
experimental measurements in predicting certain hydrocarbons and oxygenated
intermediates. Furthermore, current kinetic models are validated under
specific single pulse deposition energies, necessitating further validation
under larger deposition energy scales.

To validate the kinetics of the plasma-assisted oxidation process in $\rm
C_{3}H_{8}/O_{2}/Ar$ mixtures under large deposition energy scales, stable
discharge characteristics and spatiotemporal uniformity are required in
experiments. Fast ionization wave (FIW), a nanosecond capillary discharge at
moderate pressures, is characterized by spatiotemporal uniformity of
discharge, high reduced electric field, high electron energy, and fast
propagation speed [16]. FIW serves as an ideal platform for investigating
fundamental plasma processes under various discharge parameters and energy
deposition conditions [17, 18, 19, 20, 21, 16, 22]. Due to its uniformity and
high energy deposition, FIW is suitable for investigating plasma-assisted fuel
cracking and oxidation. FIW has two discharge types depending on the presence
or absence of a metal screen. In the absence of a metal screen (referred to as
nanosecond capillary discharge, nCD), the discharge propagation pattern shifts
from single-direction to dual-direction propagation, improving discharge
channel uniformity. Consequently, nCD is particularly suitable for validating
the accuracy of $\rm C_{3}H_{8}/O_{2}/Ar$ plasma combustion mechanisms under
varying energy deposition conditions.

In this study, we investigate the mechanism of plasma-assisted low-temperature
oxidation of propane using numerical simulations and capillary discharge
experiments. Section 2 outlines the numerical model, the kinetic schemes, and
the experimental setups. The validation of the kinetic schemes against
capillary discharge experiments and DBD experiments [13] is detailed in
Section 3. An uncertain analysis of the kinetic schemes is performed using
path analysis in Section 4.1. Section 4.2 compares the differences between nCD
and DBD in the process of propane oxidation. Finally, we provide conclusions
in Section 5.

##  2 Experiment and modelling method

###  2.1 The experiment settings

Figure 1: The shematic of the nCD assisted propane oxidation experiment.

The experimental setup of the nCD-assisted propane oxidation is depicted in
Figure 1. This setup primarily comprises a gas supply and control system, a
plasma discharge system, and a gas collection and analysis system. The
discharge system includes a high-voltage electrode, a ground electrode, a
quartz tube, and insulating rubber gaskets. The cross-sectional diameter and
length of the discharge region are 4.5 mm and 50 mm, and the distance between
the two electrodes is 44 mm.

Figure 2: The nCD experimental operation sequence diagram on the process of
gas collection.

The gas mixture in the cylinder consists of Ar, $\rm O_{2}$, and $\rm
C_{3}H_{8}$ in a ratio of 92.46:6.33:1.21. The gas flow rate is controlled
using a flowmeter (ASERT AST10-ALC 100 sccm) and is set to 30 sccm. Gas flow
is controlled using a flowmeter (ASERT AST10-ALC 100 sccm) set at 30 sccm, and
pressure is maintained at 3000 Pa by a rotary vane pump, monitored by a vacuum
gauge (AILEIKE TT&C Equipment Co., Ltd., ALKC400HJ, 0-1 atm absolute pressure
gauge, 0.05%FS). As shown in Fig. 2, valve 3 is adjusted to stabilize the gas
pressure before collecting the gas. During collection, valves 3 and 2 are
closed, and a compression tank with variable volume stabilizes and stores the
low-pressure gas, which is later pressurized to atmospheric pressure and
analyzed using a gas chromatograph (GC) (Agilent 7890B) for species
concentrations.

The discharge is powered by a nanosecond pulse discharge power supply (Xi’an
Lingfengyuan Electronic Technology Co., Ltd., HVP-20P), providing continuous
pulse outputs for steady-state gas collection. Each pulse has triangular wave
shape with 50 ns rise and fall times, operating at a frequency of 1 kHz and a
voltage range of 4.2 kV to 5 kV. Discharge voltage is monitored with a high-
voltage probe (Tektronix, P6015A), and discharge current is measured with a
current coil (Pearson 6585). Both voltage and current data are recorded using
an oscilloscope (Siglent SDS2354X 350 MHz).

###  2.2 Modelling methods

A new coupled plasma and combustion chemistry solver, CPCC, has been developed
by integrating the ZDPlasKin code [23] with a custom-built combustion code in
F90 format. CPCC supports sensitivity and path flux analysis and interfaces
with popular post-processing software such as QtPlaskin [24] and PumpKin [25].
Detailed descriptions of CPCC can be found in [26].The electron impact process
is solved using the Boltzmann equation solver BOLSIG+ [27], which is
integrated within ZDPlasKin. The plasma kinetic mechanism and ground state
kinetic mechanism are alternately solved using the splitting method. For
efficient solution of large combustion kinetic mechanisms, CPCC utilizes the
semi-implicit ODE solver VODPK [28] and the fifth-order implicit Runge-Kutta
ODE solver, RADAU5 [29].

The key parameter of the reduced electric field ($E/N$) is derived from
experimental current data using the following equation:

$I=n_{e}\mu(\frac{E}{N})EeS$ (1)

where $n_{e}$ is the electron density (initially set to $1.0\times
10^{9}cm^{-3}$), $I$ is the electron mobility dependent on the reduced
electric field, $N$ is the natural particle density of last time, $E$ is the
electric field strength, $e$ is the elementary charge, and $S$ is the cross-
sectional area of the plasma. $E$ is iteratively obtained using equation 1
through the Newton iteration method. The intial electron density is given a
empirical value, because $E$ is insensitive to the intial electron density,
which will iterate consistently to a reasonable value in the first few time
steps without affecting the result of the calculation. Then $E/N$ is obtained
by dividing $E$ by $N$.

In present model, as the Poisson equation is not solved, the deposition energy
is assumed to be an adjustable parameter and controls the discharge duration
[5]. For DBD experiment, the present model’s deposition energy stands at
$1.85\times 10^{-5}cm^{-3}$, significantly lower than the experimentally
measured value. Furthermore, the deposition energy control is not implemented.
For nCD experiment, To ensure good agreement with the experimental results,
the model utilized 40% of the experimental deposited energy to regulate the
discharge time.

To account for primary energy loss due to conduction to the quartz channel
walls, similar to the approach in Ref. [7], the model includes surface heat
transfer between the plasma and the quartz wall. In the energy equation for
the nCD experiment, an additional energy loss term is introduced:

$\rho
C_{v}\frac{\mathrm{d}T_{\mathrm{gas}}}{\mathrm{d}t}=-\sum_{i=1}^{i_{\max}}H_{i}\omega_{i}-\frac{T_{\mathrm{gas}}-T_{\mathrm{wall}}}{\left(\frac{1}{A_{1}h_{1}}+\frac{\ln\left(r_{2}/r_{1}\right)}{2\pi\lambda
L}\right)V_{p}}$ (2)

where $\rho$ is the gas density, $C_{v}$ is the specific heat capacity at
constant volume of the gas, $H_{i}$ and $\omega_{i}$ represent the enthalpy of
formation and production/loss rate of component i, respectively. $T_{wall}$ is
the outer wall temperature of the quartz tube (set to 300K), $A_{1}$ is the
contact area between the plasma and the inner wall of the quartz tube, $h_{1}$
is the heat transfer coefficient between the plasma and the inner wall of the
quartz tube (which has been adjusted to 10 $W/(m^{2}K)$ based on experimental
results due to the cylindrical shape and small diameter of this structure,
making empirical formulas unsuitable for estimating the heat transfer
coefficient), $r_{2}$ and $r_{1}$ are the outer and inner radius of the quartz
tube, $\lambda$ is the thermal conductivity of quartz (which is 1.09
$W/(mK)$), $L$ is the length of the discharge region, and $V_{p}$ is the
volume of the plasma.

For a DBD configuration where the discharge region is rectangular, the energy
loss term in the gas energy equation is adjusted as follows:

$\rho
C_{v}\frac{\mathrm{d}T_{\mathrm{gas}}}{\mathrm{d}t}=-\sum_{i=1}^{i_{\max}}H_{i}\omega_{i}-\frac{T_{\mathrm{gas}}-T_{\mathrm{wall}}}{\left(\frac{1}{A_{2}h_{2}}+\frac{\delta}{A_{2}\lambda}\right)V_{p}}$
(3) $h_{2}=\frac{\lambda}{LN_{u}}$ (4) $N_{u}=0.664R_{e}^{1/2}P_{r}^{1/3}$ (5)
$P_{r}=\frac{C_{p}\mu}{\lambda}$ (6) $R_{e}=\frac{\rho vL}{\mu}$ (7)

where $A_{2}$ is the contact area between the plasma and the upper and lower
walls, $h_{2}$ is the surface heat transfer coefficient of the walls
(calculated using equation 4), and the Nusselt number $N_{u}$ is calculated
using equation 5 from the book [30]. $P_{r}$ is the Prandtl number (calculated
by equation 6), $R_{e}$ is the Reynolds number (calculated by equation 7),
$\mu$ is the dynamic viscosity, $v$ is the airflow velocity, and $C_{p}$ is
the specific heat capacity at constant pressure of the gas.

This kinetic model includes both plasma chemistry and classical combustion
chemistry. The combustion chemistry model incorporates a total of 2521
reactions, which have been updated based on the scheme proposed in Ref. [13].
Specific modifications derived from experimental data are outlined in Table 1.
In addition, several combustion reactions have been added and replaced,
drawing from various established mechanisms such as HP-Mech v3.3 [31],
AramoMech3.0 [32], NUIGMech1.1 [33] (see supplementary material of Table S1).
The plasma chemistry model consists of 2202 reactions, involving vibrationally
excited species such as $\rm C_{3}H_{8}(v1)$-$\rm C_{3}H_{8}(v27)$, $\rm
CH_{4}(v2,4)$, $\rm CH_{4}(v1,3)$, $\rm O_{2}(v1)$-$\rm O_{2}(v4)$;
electronically excited species like $\rm O_{2}(a^{1}\Delta_{g})$, $\rm
O_{2}(b^{1}\Sigma_{g}^{+})$, $\rm O_{2}(4.5eV)$(summation of $\rm
O_{2}\left(c^{1}\Sigma_{u}^{-}\right)$, $\rm
O_{2}\left(C^{3}\Delta_{u}\right)$, $\rm
O_{2}\left(A^{3}\Sigma_{u}^{+}\right)$), $\rm O(^{1}D)$, $\rm O(^{1}S)$, $\rm
Ar(1s_{2})$, $\rm Ar(1s_{4})$, $\rm Ar(11.55eV)$, $\rm Ar(2p_{1})$; ions like
$\rm O^{+}$, $\rm O_{2}^{+}$, $\rm O^{-}$, $\rm O_{2}^{-}$, $\rm O_{3}^{-}$,
$\rm CH_{4}^{+}$, $\rm CH_{3}^{+}$, $\rm CH_{2}^{+}$, $\rm CH_{5}O^{+}$, $\rm
C_{2}H_{3}O^{+}$, $\rm CHO_{3}^{-}$, $\rm CO_{3}^{-}$, $\rm CHO_{2}^{-}$, $\rm
OH^{-}$, $\rm C_{2}H_{4}^{+}$, $\rm C_{2}H_{5}^{+}$, $\rm C_{2}H^{+}$, $\rm
C_{2}H_{3}^{+}$, $\rm H_{5}C_{2}O^{+}$, $\rm C_{3}H_{7}^{+}$, $\rm
C_{4}H_{9}^{+}$, $\rm C_{3}H_{8}^{+}$, and electron. The plasma chemistry
model includes chemical reactions involving vibrationally excited states,
vibrational relaxation reactions, charge transfer reactions, and reaction
channels involving $\rm O(^{1}D)$, $\rm O(^{1}S)$, $\rm
O_{2}(a^{1}\Delta_{g})$ , $\rm O_{2}(b^{1}\Sigma_{g}^{+})$, $\rm
O_{2}(4.5eV)$, $\rm Ar(1s_{2})$, $\rm Ar(1s_{4})$, $\rm Ar(11.55eV)$, $\rm
Ar(2p_{1})$ with fuel fragments. The electron collision dissociation cross-
section for $\rm C_{3}H_{8}$ is sourced from Ref. [34], while reactions
involving electron collision dissociation for $\rm CH_{4}$, including
dissociation and ionization, are adopted from Ref. [6]. The plasma chemistry
of Ar and $\rm O_{2}$ are taken from Ref. [35].

Table 1: Updated plasma reactions and combustion reactions in this work.

| Plasma reaction |  | Original rate
---
constant
| Updated rate
---
constant
R1 | $\rm e+C_{3}H_{8}\rightarrow e+H_{2}+C_{3}H_{6}$ |  $f$($\sigma_{1}$)[34] |  $0.8\times f$($\sigma_{1}$)
R2 | $\rm e+C_{3}H_{8}\rightarrow e+CH_{2}+C_{2}H_{6}$ |  $f$($\sigma_{2}$)[34] |  $12\times f$($\sigma_{2}$)
R3 | $\rm e+C_{3}H_{8}\rightarrow e+CH_{4}+C_{2}H_{4}$ |  $f$($\sigma_{3}$)[34] |  $0.36\times f$($\sigma_{3}$)
R4 | $\rm e+C_{3}H_{8}\rightarrow e+H_{2}+H_{2}+C_{3}H_{4}$ |  $f$($\sigma_{4}$)[34] |  $0.75\times f$($\sigma_{4}$)
R5 | $\rm e+C_{3}H_{8}\rightarrow e+CH_{2}+C_{2}H_{2}+2H_{2}$ | 0.0 |  $21\times f$($\sigma_{2}$)
| Combustion reaction |  | Original Arrhenius
---
equation of A
| Updated Arrhenius
---
equation of A
R6 | $\rm CH_{3}OH+O_{2}\leftrightarrow CH_{3}O+HO_{2}$ | $3.58\times 10^{4}$ |  0.0 [36]
R7 |  $\rm n$-$\rm C_{3}H_{7}O_{2}+HO_{2}\leftrightarrow n$-$\rm C_{3}H_{7}O_{2}H+O_{2}$ | $1.75\times 10^{10}$ | $6.0\times 10^{8}$
R8 |  $\rm i$-$\rm C_{3}H_{7}O_{2}+HO_{2}\leftrightarrow i$-$\rm C_{3}H_{7}O_{2}H+O_{2}$ | $1.75\times 10^{10}$ | $6.0\times 10^{8}$
R9 |  $\rm i$-$\rm C_{3}H_{7}O_{2}+C_{2}H_{5}O_{2}\leftrightarrow CH_{3}COCH_{3}+C_{2}H_{5}OH+O_{2}$ | 0.0 | $5.0\times 10^{16}$
R10 |  $\rm n$-$\rm C_{3}H_{7}O_{2}+n$-$\rm C_{3}H_{7}O_{2}\leftrightarrow C_{2}H_{5}CHO+n$-$\rm C_{3}H_{7}OH+O_{2}$ | $1.42\times 10^{16}$ | $7.0\times 10^{15}$
R11 |  $\rm n$-$\rm C_{3}H_{7}O_{2}+i$-$\rm C_{3}H_{7}O_{2}\leftrightarrow C_{2}H_{5}CHO+i$-$\rm C_{3}H_{7}OH+O_{2}$ | $1.409\times 10^{16}$ | $4.109\times 10^{16}$
R12 |  $\rm i$-$\rm C_{3}H_{7}O_{2}+n$-$\rm C_{3}H_{7}O_{2}\leftrightarrow CH_{3}COCH_{3}+n$-$\rm C_{3}H_{7}OH+O_{2}$ | $1.409\times 10^{16}$ | $1.169\times 10^{15}$
R13 |  $\rm i$-$\rm C_{3}H_{7}O_{2}+i$-$\rm C_{3}H_{7}O_{2}\leftrightarrow CH_{3}COCH_{3}+i$-$\rm C_{3}H_{7}OH+O_{2}$ | $1.409\times 10^{16}$ | $2.109\times 10^{16}$

  * *

$\sigma_{x}$ indicate the cross sections for the electron-impact reaction. The
rate constant values are calculated by BOLSIG+ (indicated by $f$ ) given E/N
and the cross section data. The unit of cross section data is $m^{2}$. The
unit of arrhenius equation pre-exponential factor A is
$cm^{3}\>mol^{-1}\>s^{-1}$.

The original model in Section 3 and Section 4 is defined as the plasma
chemistry and classical combustion chemistry in which the reaction rate
coefficients are from literature without any modifications. Compared to the
model of Ref. [13], the original model reconstruct the plasma chemistry and
has the same combustion chemistry. As some reactions are not validated with
experiments or quantum chemical calculations, the original model predictions
for some species have large uncertainties. So the electron impact cross
sections of propane, branch ratios of self-/cross-reactions for $\rm i$-/$\rm
n$-$\rm C_{3}H_{7}O_{2}$ and H-abstraction of $\rm i$-/$\rm n$-$\rm
C_{3}H_{7}O_{2}$ are fitted within uncertainties according to experimental
data and path flux analysis (named as updated model in Section 3 and Section
4). The updated reactions are listed in Table 1. Detailed descriptions and
uncertainty analysis can be found in Section 4.1. The species overall
predictions of the updated model are highly improved compared to the original
model in DBD experiment [13] and nCD experiment conditions. Even so, the
quantum chemical calculations or fundamental kinetic experiment are still
requred to improve the model’s prediction capacity.

##  3 Validation of the mechanism

###  3.1 Validation I:DBD assisted propane oxidation

Figure 3: Comparison between measurement and predictions from different
kinetics schemes species in repetitively-pulsed nanosecond DBD. The deep blue
and light blue bars represent experiment measurement and model prediction
results from the Ref. [13].

Model validation is performed by comparing the model predicted and
experimental molar fractions of speices in a repetitively-pulsed nanosecond
DBD [13], as shown in Fig. 3. The experimental conditions involve a deposition
energy of $1.64\times 10^{-4}J/cm^{3}$, at 30 Torr and a discharge region
temperature of 340 K, with a $\rm C_{3}H_{8}/O_{2}/Ar$ mixture (mole fraction
4:20:76). The updated model’s predictions for several species, including $\rm
H_{2}O$, CO, $\rm CH_{2}O$, $\rm CH_{3}OH$, $\rm C_{2}H_{2}$, $\rm
C_{2}H_{4}$, $\rm C_{2}H_{5}OH$, $\rm C_{2}H_{5}OOH$, $\rm C_{3}H_{4}$-$\rm
A$, $\rm C_{3}H_{4}$-$\rm P$, $\rm C_{3}H_{6}$, $\rm CH_{3}COCH_{3}$, $\rm
C_{2}H_{5}CHO$, $\rm i$-$\rm C_{3}H_{7}OH$, and $\rm C_{3}H_{7}OOH$ are found
to be in good agreement with experimental data, falling within the range of
experimental measurement errors. The updated model’s predictions for species
$\rm CO_{2}$, $\rm CH_{4}$, $\rm C_{2}H_{6}$, $\rm CH_{2}CO$, $\rm CH_{3}CHO$,
and $\rm n$-$\rm C_{3}H_{7}OH$ are within 5 times of the experimental
measurement. However, the updated model’s predictions for species $\rm
CH_{3}$, $\rm CH_{3}OOH$, $\rm cC_{3}H_{6}$, $\rm C_{2}H_{3}CHO$, $\rm
C_{3}H_{6}O_{1-3}$, and $\rm C_{3}H_{6}O_{1-2}$ are over an order of magnitude
off. The updated model demonstrates noticeable improvements in prediction
accuracy for many species compared to the original model and a previously
published model [13].

The reasons for the lower predictions of $\rm C_{3}H_{6}O_{1\text{-}2}$ and
$\rm C_{3}H_{6}O_{1\text{-}3}$ in the updated model are as follows: In Ref.
[37, 38], $\rm C_{3}H_{6}$ is consumed by an O atom to form $\rm
C_{3}H_{6}O_{1\text{-}2}$ and $\rm C_{3}H_{6}O_{1\text{-}3}$ intermediates
through the addition of an O atom to the C=C bonds in the center and end
carbon atoms of $\rm C_{3}H_{6}$. Then, $\rm C_{3}H_{6}O_{1\text{-}2}$ and
$\rm C_{3}H_{6}O_{1\text{-}3}$ intermediates decompose into $\rm CH_{3}$, $\rm
C_{2}H_{5}$, $\rm CH_{2}O$, and H, among others. The reactions provided in
Ref. [37] do not consider these intermediate products ($\rm
C_{3}H_{6}O_{1\text{-}2}$ and $\rm C_{3}H_{6}O_{1\text{-}3}$) but instead
directly provide the reactions leading to the final products. This ultimately
results in the underprediction of $\rm C_{3}H_{6}O_{1\text{-}2}$ and $\rm
C_{3}H_{6}O_{1\text{-}3}$.

###  3.2 Validation II: nCD assisted propane oxidation

The DBD experiment was conducted under conditions with a deposition energy of
$1.64\times 10^{-4}J/cm^{3}$. To further validate the accuracy of the
mechanism, nCD experiments were conducted with deposition energies ranging
from $1.7\times 10^{-4}J/cm^{3}$ to $1.53\times 10^{-3}J/cm^{3}$. The
uniformity of the nCD was evaluated using images captured by a high-speed ICCD
camera operated at 5.0 kV, which were processed using an Abel transformator,
as depicted in Fig. 4. Due to the enhanced electric field near the electrode
tip, localized discharge enhancement was observed on the ground electrode.
Conversely, a phenomenon of discharge weakening was noted at distances ranging
from 9 to 20 mm away from the high-voltage electrode. These localized effects
were considered negligible in comparison to the overall discharge
characteristics. Overall, the nCD was found to exhibit relatively uniform
behavior across the tested conditions.

Figure 4: The nCD image captured by high-speed ICCD camera and abel
transformed image with voltage of 5.0 kV. The gate width is set to 100ns.

Figure 5 presents a comparison between experimental measurements and
predictions from different kinetic schemes in the nCD conditions. The original
model exhibits poor predictions for species such as CO, $\rm H_{2}$, $\rm
CO_{2}$, $\rm CH_{4}$, $\rm C_{2}H_{6}$, $\rm C_{2}H_{2}$, $\rm C_{3}H_{6}$.
In contrast, the updated model shows significant improvements across these
species. However, despite these enhancements, notable discrepancies remain in
the predicted mole fractions of $\rm CO_{2}$, $\rm C_{2}H_{4}$, and CO, which
are under-predicted by 55%, 45%, and 45%, respectively.

Figure 5: Comparison between measurement and predictions from different
kinetics scheme in nCD at different voltages.

The comparison between measurement and predictions from different kinetics
scheme in nCD is shown in Fig. 5. The original model has poor predictions in
species CO, $\rm H_{2}$, $\rm CO_{2}$, $\rm CH_{4}$, $\rm C_{2}H_{6}$, $\rm
C_{2}H_{2}$, $\rm C_{3}H_{6}$. The updated model has great improvement in
above mentioned species. However, compared to experimental measurements, the
updated model still has significant differences in the mole fractions of $\rm
CO_{2}$, $\rm C_{2}H_{4}$, and CO. The mole fractions of updated model are
under-predicted by 55%, 45%, and 45%, respectively. The predictions of
remaining species are generally well.

The updated model’s predictions for $\rm C_{3}H_{8}$, $\rm H_{2}$, $\rm
CH_{4}$, and $\rm C_{2}H_{6}$ generally fall within the experimental error
range. For $\rm C_{3}H_{6}$, while it is overpredicted by 30% at low
deposition energy, the agreement improves at higher deposition energies. The
trend of $\rm C_{2}H_{2}$ predicted by the updated model aligns with
experimental data, albeit with molar fractions approximately 30% lower than
measured values. This suggests that while the added reaction pathways in the
updated model are correct, further research is needed to refine the cross-
sectional data used.

##  4 Results and discussion

###  4.1 Uncertain analysis in the plasma and combustion kinetic mechanism

The updated model was derived from the original model through adjustments to
several reactions based on experimental findings and pathway flux analysis.
Specifically, the modifications listed in Table 1 are discussed in detail as
follows.

Figure 6: The pathway analysis of $\rm i\text{-}/n\text{-}C_{3}H_{7}OOH$ in
the updated model and the original model on the process of DBD assisted
propane oxidation.

In the original model, significant over-predictions of $\rm
i\text{-}/n\text{-}C_{3}H_{7}OOH$ were observed during the DBD-assisted
propane oxidation process, as depicted in Fig. 3. The generation rate of $\rm
i\text{-}/n\text{-}C_{3}H_{7}OOH$ was found to be two orders of magnitude
larger than its consumption rate. Pathway flux analysis based on the original
model identified that the primary generation pathways for $\rm
i\text{-}/n\text{-}C_{3}H_{7}OOH$ involved reactions $\rm
i\text{-}/n\text{-}C_{3}H_{7}O_{2}$ with $\rm HO_{2}$ (reactions R7 and R8,
respectively), as illustrated in Fig. 6. However, the rate coefficients for R7
and R8 used in the original model were adapted from analogous reactions such
as those between $\rm CH_{3}O_{2}$ and $\rm HO_{2}$ reported in Ref. [39].
These rate coefficients are valid within a temperature range below 338 K and
are subject to an error of up to 5. Given the temperature range of 340-500 K
in the current study and the inherent uncertainties associated with the
analogy method, significant uncertainties were anticipated for R7 and R8.
Therefore, adjustments were made to the rate coefficients of R7 and R8,
reducing them by a factor of 0.034 based on insights from the DBD experiment
[13]. The updated model now indicates that 87% of $\rm n\text{-}C_{3}H_{7}OOH$
is generated from the reaction of $\rm n\text{-}C_{3}H_{7}O$ with OH, while
84% of $\rm i\text{-}C_{3}H_{7}OOH$ is produced via reaction R8.

Figure 7: The pathway analysis of $\rm CH_{3}OH$, $\rm CH_{3}O$ and $\rm
CH_{2}O$ in the updated model and the original model on the process of DBD
assisted propane oxidation.

Both the model by Chen et al [13] and the original model exhibit under-
prediction of $\rm CH_{2}O$ and over-prediction of $\rm CH_{3}OH$. As depicted
in Fig. 7, pathway flux analysis of the original model indicates that 97% of
$\rm CH_{3}OH$ is generated via the reaction of $\rm CH_{3}O$ with $\rm
HO_{2}$ (R6). However, it was noted that the reaction rate coefficients of R6
were sourced from Ref. [36], which states that the reaction between $\rm
CH_{3}OH$ and $\rm O_{2}$ produces $\rm CH_{2}OH$ and $\rm HO_{2}$, not $\rm
CH_{3}O$ and $\rm HO_{2}$. Consequently, in the updated model, reaction R6 has
been removed. The updated model’s prediction for $\rm CH_{3}OH$ now aligns
well with experimental measurements, with 63% of $\rm CH_{3}OH$ being produced
by the self-reaction of $\rm CH_{3}O$ (Reaction 8), and 34% being generated
through Reaction 9.

$\rm CH_{3}O+CH_{3}O\leftrightarrow CH_{3}OH+CH_{2}O$ (8) $\rm
CH_{3}O_{2}+OH\leftrightarrow CH_{3}OH+O_{2}$ (9)

To address the lack of certain $\rm CH_{3}O$ consumption pathways, the $\rm
CH_{3}O$ sub-mechanism was supplemented by referencing the HP-Mech v3.3
mechanism [31] (as shown in Table LABEL:Combustion_reactions). The removal of
R6 and the supplement of $\rm CH_{3}O$ consumption reactions led to
significant changes in the $\rm CH_{3}O$ consumption pathway, as illustrated
in Fig. 7. In the original model, $\rm CH_{3}O$ was almost exclusively
consumed by the reverse reaction of R6. In contrast, in the updated model,
half of $\rm CH_{3}O$ is consumed to form $\rm CH_{3}$, while the remaining
$\rm CH_{3}O$ is consumed to form $\rm CH_{2}O$ and $\rm CH_{2}O$.
Consequently, the prediction of $\rm CH_{2}O$ in the updated model has
improved and now aligns well with experimental measurements due to the removal
of R6 and the supplement of $\rm CH_{3}O$ consumption reactions.

$\rm n\text{-}C_{3}H_{7}O\rightarrow C_{2}H_{5}+CH_{2}O$ (10)
$\mathrm{C}_{3}\mathrm{H}_{5}\text{-}\mathrm{T}+\mathrm{O}_{2}\leftrightarrow\mathrm{CH}_{3}\mathrm{CO}+\mathrm{CH}_{2}\mathrm{O}$
(11)
$\mathrm{CH}_{3}+\mathrm{O}\leftrightarrow\mathrm{CH}_{2}\mathrm{O}+\mathrm{H}$
(12) $\rm CH_{3}O+O_{2}/H/OH/O\leftrightarrow CH_{2}O+HO_{2}/H_{2}/H_{2}O/OH$
(13)

Figure 8: The pathway analysis of $\rm n\text{-}C_{3}H_{7}OH$ and $\rm
i\text{-}C_{3}H_{7}OH$ in the updated model and the original model on the
process of DBD assisted propane oxidation.

As depicted in Fig. 3, the original model exhibited over-predictions of $\rm
n\text{-}C_{3}H_{7}OH$ and under-predictions of $\rm i\text{-}C_{3}H_{7}OH$.
It was observed that nearly all $\rm i\text{-}C_{3}H_{7}OH$ and almost half of
$\rm n\text{-}C_{3}H_{7}OH$ were produced through self- and cross-reactions of
$\rm n\text{-}/i\text{-}C_{3}H_{7}O_{2}$ (Reactions R10-R13), as illustrated
in Fig. 8. The inaccuracies in the branching ratios of these reactions likely
contributed to the over-prediction of $\rm n\text{-}C_{3}H_{7}OH$ and the
under-prediction of $\rm i\text{-}C_{3}H_{7}OH$. The rate coefficients for
Reactions R10-R13 were originally derived from Ref. [9]. However, the
temperature-dependent rate constants used in Ref. [9], sourced from Ref. [40],
lacked detailed pressure and temperature dependencies. Moreover, the products
of Reactions R10-R13 were not experimentally validated in Ref. [9],
contributing to significant uncertainty in the branch ratios. To improve
prediction accuracy, this model adjusted the branching ratios of Reactions
R10-R13 based on DBD experimental results. This adjustment notably enhanced
the prediction accuracy for $\rm i\text{-}C_{3}H_{7}OH$, although a slight
increase in $\rm n\text{-}C_{3}H_{7}OH$ mole fraction compared to the original
model was observed. This increase can be attributed to the enhanced reaction
rate of $\rm n\text{-}C_{3}H_{7}O$ with $\rm HO_{2}$, reflecting the increased
mole fraction of $\rm HO_{2}$ resulting from the deletion of Reaction R6.

Additionally, a formation pathway for $\rm C_{2}H_{5}OH$ was introduced based
on the NUIGMech1.1 mechanism [33], using rate coefficients derived from self-
and cross-reactions similar to those of $\rm
n\text{-}/i\text{-}C_{3}H_{7}O_{2}$. This inclusion brought the model’s
prediction for $\rm C_{2}H_{5}OH$ within the experimental measurement error
range, as depicted in Fig. 3.

Furthermore, the original model showed under-prediction of $\rm C_{2}H_{5}OOH$
compared to experimental measurements, with almost all $\rm C_{2}H_{5}OOH$
formed via the reaction between $\rm C_{2}H_{5}O_{2}$ and $\rm HO_{2}$. The
increased mole fraction of $\rm HO_{2}$ in the updated model led to an
enhanced reaction rate between $\rm C_{2}H_{5}O_{2}$ and $\rm HO_{2}$,
aligning the predicted mole fraction of $\rm C_{2}H_{5}OOH$ more closely with
experimental measurements. This further supports the decision to delete
Reaction R6.

Figure 9: The pathway analysis of $\rm CH_{4}$, $\rm C_{2}H_{2}$, $\rm
C_{2}H_{4}$, $\rm C_{2}H_{6}$ and $\rm C_{3}H_{6}$ in the updated model and
the original model on the process of nCD assisted propane oxidation.

The discrepancies of $\rm CH_{4}$, $\rm C_{2}H_{2}$, $\rm C_{2}H_{4}$, $\rm
C_{2}H_{6}$ and $\rm C_{3}H_{6}$ between experimental measurements and the
original model are observed in the DBD experiment and the nCD experiment.
Given the large energy deposition scales of the nCD experiment, it was chosen
as the primary reference for adjusting the kinetic model. A discharge peak
voltage of 5.0 kV was selected for pathway flux analysis, as shown in Fig. 9.
In the original model, $\rm CH_{4}$, $\rm C_{2}H_{4}$, $\rm C_{2}H_{6}$, and
$\rm C_{3}H_{6}$ are produced from the electron impact dissociation of $\rm
C_{3}H_{8}$. The dissociation cross sections for $\rm C_{3}H_{8}$ were
obtained from Ref. [34], which provided estimates based on the total
experimental dissociation cross sections for both ionized and neutral
fragments ($\rm e+C_{3}H_{8}\rightarrow A+B+e$ and $\rm
e+C_{3}H_{8}\rightarrow A^{+}+B+2e$). However, the cross-section branch ratios
were estimated using empirical formulas, leading to significant uncertainties.
Consequently, the reaction rate coefficients for $\rm C_{3}H_{8}$ dissociation
were adjusted according to experimental results.

First, in the original model, Reaction (R1) is the primary production pathway
for $\rm C_{3}H_{6}$, and its prediction was overestimated across all
discharge voltage scales in the nCD experiment. To improve accuracy, the cross
section for Reaction (R1) was reduced by a factor of 0.8, which enhanced the
prediction accuracy for $\rm C_{3}H_{6}$, as shown in Fig. 5.

Second, reaction (R2) is responsible for the formation of $\rm C_{2}H_{6}$.
The original model significantly under-predicted $\rm C_{2}H_{6}$ in the nCD
experiment. An adjustment factor of 12 was applied to the cross section for
Reaction (R2), aligning the $\rm C_{2}H_{6}$ predictions with the experimental
results.

Furthermore, reaction (R3) is the primary pathway for the formation of both
$\rm CH_{4}$ and $\rm C_{2}H_{4}$. The plasma sub-mechanism for $\rm CH_{4}$
was adopted from Ref. [6, 5] with experimental validation, indicating that the
$\rm CH_{4}$ consumption rates in the original model were accurate. Despite
this, both $\rm CH_{4}$ and $\rm C_{2}H_{4}$ were over-predicted in the DBD
experiment, and $\rm CH_{4}$ was also over-predicted in the nCD experiment.
This suggested that the cross section for Reaction (R3) was overestimated. An
adjustment to the cross section for Reaction (R3) was made to achieve better
agreement for $\rm CH_{4}$ in the nCD experiment and $\rm C_{2}H_{4}$ in the
DBD experiment.

Besides, A discrepancy in the prediction of $\rm C_{3}H_{4}$ between the
original model and the DBD experiment was noted. Reaction (R4) is responsible
for 93% of the formation of $\rm C_{3}H_{4}\text{-}A$ and $\rm
C_{3}H_{4}\text{-}P$. Thus, the rate of Reaction (R4) was adjusted to better
match the DBD experimental results.

Finally, the original model under-predicted $\rm C_{2}H_{2}$ in both the nCD
and DBD experiments, as shown in Fig. 5 and Fig. 3. In the original model,
$\rm C_{2}H_{2}$ is generated through reactions between excited argon and $\rm
C_{2}H_{4}$ and between $\rm O$ and $\rm C_{3}H_{4}\text{-}A$. However, with a
mole fraction of $\rm C_{2}H_{2}$ around 600 ppm in the nCD experiment
compared to a $\rm C_{2}H_{4}$ mole fraction of less than 400 ppm, $\rm
C_{2}H_{4}$ alone cannot account for the $\rm C_{2}H_{2}$ formation.
Similarly, $\rm C_{3}H_{4}\text{-}A$ is present at too low a concentration to
significantly contribute to $\rm C_{2}H_{2}$ production. Therefore, electron-
impact dissociation of $\rm C_{3}H_{8}$ (Reaction R5) was considered as a
primary source for $\rm C_{2}H_{2}$. Adding Reaction (R5) to the updated model
significantly improved the accuracy of $\rm C_{2}H_{2}$ predictions, with the
updated model’s $\rm C_{2}H_{2}$ predictions closely matching experimental
results from both the nCD and DBD experiments.

As stated above, The model undergoes significant updates through pathway flux
analysis informed by experiments conducted under nCD and DBD conditions. By
adjusting key reaction mechanisms according to experimental data, including
$\rm C_{3}H_{6}$, $\rm C_{2}H_{6}$, $\rm CH_{4}$, $\rm C_{2}H_{4}$, $\rm
C_{2}H_{2}$, and $\rm C_{3}H_{4}$, the updated model now demonstrates markedly
improved predictive accuracy. This enhancement is particularly evident in its
ability to better replicate experimental observations in both nCD and DBD
experiments.

###  4.2 The comparisons of $\rm C_{3}H_{8}$ oxidation process assisted by DBD
and nCD

To explore the impact of the single pulse deposition energy (SPDE) on the
oxidation of $\rm C_{3}H_{8}$, we have designed three distinct cases, each
tailored to mimic specific energy deposition conditions. As depicted in Table
2, these cases are denoted as Case1, Case2, and Case3, with corresponding SPDE
values in the model set at $1.85\times 10^{-5}$ $J/cm^{3}$, $9.08\times
10^{-5}$ $J/cm^{3}$, and $3.9\times 10^{-4}$ $J/cm^{3}$, respectively.
Notably, the total deposition energy in the model is standardized at
$7.36\times 10^{-3}$ $J/cm^{3}$ by adjusting the discharge frequency and pulse
number, ensuring consistency across the cases. Furthermore, the gas
temperature is maintained at a constant 400 K, while the discharge current
values for the DBD and nCD in the models are sourced from Ref. [13] and the
nCD experiment, respectively. These cases are designed to provide insights
into the role of SPDE on $\rm C_{3}H_{8}$ oxidation, offering valuable
comparisons under varying energy deposition conditions.

Table 2: Initial parameters of simulation cases.

| Case1 | Case2 | Case3
---|---|---|---
Discharge type | DBD | nCD | nCD
Peak voltage, kV | 9.0 | 4.5 | 5.0
Discharge time, ms | 20
Frequency, Hz | 20000 | 4050 | 950
Pulse number | 400 | 81 | 19
Single pulse deposition energy in model, $J/cm^{3}$ | $1.85\times 10^{-5}$ | $9.08\times 10^{-5}$ | $3.9\times 10^{-4}$
Total deposition energy in model, $J/cm^{3}$ | $7.36\times 10^{-3}$
Gas temperature, K | 400
Pressure, Pa | 4000
Discharge mixture |  $\rm C_{3}H_{8}/O_{2}/Ar$ 4:20:76

In Fig. 10(a), the average reduced electric field evolution for DBD and nCD is
depicted, where the reduced electric field is averaged over all discharge
pulses. Notably, the average reduced electric field for nCD experiences a
rapid increase to 450-500 Td before 3 ns, followed by a dramatic decrease to
150-200 Td at 8 ns, maintaining this level until the end of the discharge. In
contrast, the average reduced electric field for DBD fluctuates repeatedly
from 0 Td to 400 Td. As illustrated in Fig. 10(b), the reduced electric field
distribution for nCD is concentrated at 150-200 Td and 450-500 Td. In
comparison, the reduced electric field distribution for DBD predominantly
ranges from 0-25 Td and 50-250 Td, accounting for approximately 55% and 35%,
respectively, indicating a lower reduced electric field distribution compared
to nCD.

Figure 10: The comparison of average reduced electric field and reduced
electric field distribution in DBD and nCD.

Fig. 11 presents a comparison of species density in the $\rm C_{3}H_{8}$
oxidation process assisted by DBD and nCD. Interestingly, the density of
species such as $\rm H_{2}O$, CO, $\rm CO_{2}$, $\rm C_{2}H_{2}$, and $\rm
C_{2}H_{4}$ in the $\rm C_{3}H_{8}$ oxidation process remains similar between
DBD and nCD, despite varying SPDE values. This suggests that SPDE has minimal
influence on the $\rm C_{3}H_{8}$ oxidation process when the total deposition
energy remains constant. Furthermore, Fig. 12 highlights the time evolution of
$\rm C_{3}H_{6}$ and $\rm CH_{4}$ in the $\rm C_{3}H_{8}$ oxidation process
assisted by DBD and nCD. Both $\rm C_{3}H_{6}$ and $\rm CH_{4}$ are produced
during nanosecond pulse discharge and subsequently consumed in the discharge
gaps. Notably, the increase in $\rm C_{3}H_{6}$ and $\rm CH_{4}$ during a
given period correlates with the rise in SPDE. These observations collectively
suggest that the $\rm C_{3}H_{8}$ oxidation process is primarily controlled by
total deposition energy and displays limited correlation with the specific
discharge type, be it DBD or nCD.

Figure 11: Comparison of species denisty in $\rm C_{3}H_{8}$ oxidation process
assisted by DBD and nCD.

Figure 12: $\rm C_{3}H_{6}$ and $\rm CH_{4}$ time evolutions in $\rm
C_{3}H_{8}$ oxidation process assisted by DBD and nCD.

As depicted in Fig. 13, the consumption pathways of $\rm C_{3}H_{8}$ in DBD
and nCD environments were analyzed. In DBD, 57% of $\rm C_{3}H_{8}$ is
consumed through vibrationally excited reactions, attributed to the
concentrated reduced electric field distribution of 0-25 Td, which enhances
efficient vibrational excitation. Conversely, in nCD, a lower proportion of
$\rm C_{3}H_{8}$ is consumed via vibrationally excited reactions, with a
higher contribution from species such as $\rm O(^{1}D)$, $\rm Ar^{+}$, $\rm
e$, and OH. $\rm C_{3}H_{8}(v)$ is quenched back to $\rm C_{3}H_{8}$ through
relaxation reactions, contributing thermally to the process. Although a
significant portion of $\rm C_{3}H_{8}$ in DBD is excited to $\rm
C_{3}H_{8}(v)$, the energy deposited into $\rm C_{3}H_{8}(v)$ remains
relatively low due to the limited vibrational excitation energy. In DBD, the
majority of energy deposition occurs within the 100-250 Td range, similar to
processes occurring at 150-200 Td and 500 Td in nCD. Consequently, the species
density in the $\rm C_{3}H_{8}$ oxidation process, whether assisted by DBD or
nCD, exhibits comparable values. Additionally, the time evolutions of $\rm
C_{3}H_{6}$ and $\rm CH_{4}$ densities in both nCD and DBD environments
display similar tendencies.

Figure 13: $\rm C_{3}H_{8}$ consumption pathways analysis in DBD and nCD.

##  5 Conclusions

This study investigates the kinetics of plasma-assisted propane oxidation in
nanosecond pulse discharge. A zero-dimensional plasma-combustion coupled
solver equipped with sensitivity and pathway flux analysis modules was
developed, based on the ZDPlaskin package and a custom-built combustion
chemistry code. A kinetic model was proposed and validated against the nCD
experiment using gas chromatography and the nanosecond DBD experiment using
synchrotron photoionization mass spectrometry. The updated model demonstrates
predictions for 15 species falling within the range of DBD experimental
measurement errors, with significantly improved prediction accuracy for 17
species. Notably, the updated model exhibits good prediction accuracy with the
nCD experiment for $\rm C_{3}H_{8}$, $\rm H_{2}$, $\rm CH_{4}$, and $\rm
C_{2}H_{6}$ under large deposition energy scales. However, the model-predicted
fractions of $\rm CO_{2}$, $\rm C_{2}H_{4}$, and CO are 55%, 45%, and 45%
lower than those measured in the nCD experiment, respectively.

Several modifications to the reactions were made based on experiments and
pathway flux analysis. This included fitting the reaction rate coefficients of
$\rm i\text{-}/n\text{-}C_{3}H_{7}O_{2}$ with $\rm HO_{2}$ (R7, R8) due to the
large uncertainty of the analogy method. Additionally, the deletion of the
erroneous addition reaction $\rm CH_{3}OH+O_{2}\leftrightarrow CH_{3}O+HO_{2}$
resulted in improved model predictions for $\rm CH_{3}OH$. The supplementation
of the $\rm CH_{3}O$ sub-mechanism from the HP-Mech v3.3 mechanism enhanced
the updated model’s prediction accuracy for $\rm CH_{2}O$. Furthermore,
adjustments to the $\rm C_{3}H_{8}$ dissociation reaction rate coefficients,
owing to large uncertainties, greatly improved the prediction accuracy for
$\rm H_{2}$, $\rm CH_{4}$, $\rm C_{2}H_{6}$, $\rm C_{2}H_{2}$, and $\rm
C_{3}H_{6}$.

Comparisons of $\rm C_{3}H_{8}$ oxidation processes assisted by DBD and nCD
under varying SPDE but with the same total deposition energy revealed that the
process is primarily controlled by total deposition energy, demonstrating
minimal influence from SPDE and displaying limited correlation with the
specific discharge types, DBD and nCD.

Overall, these findings illustrated the intricate interplay between plasma and
combustion chemistry, offering valuable insights into the mechanisms of
plasma-assisted propane oxidation.

## Acknowledgments

The work is supported by the National Natural Science Foundation of China
(No.52277168, No.52025064, No.52207151), the National Key Research and
Development Program of China (No.2023YFB4005700, No.2023YFB4005705, and
No.2023YFB4005702-03), the University Synergy Innovation Program of Anhui
Province (No.GXXT-2022025), and the independent project of the Energy Research
Institute of Hefei Comprehensive National Science Center (Anhui Energy
Laboratory) (No.22KZZ525, No.23KZS402, No.22KZS301, and No.22KZS304). The
authors are thankful to the young research group in Atelier des Plasmas for
fruitful discussions.

## References

  * Ju and Sun [2015] Y. Ju, W. Sun,  Plasma assisted combustion: Dynamics and chemistry,  PROG ENERG COMBUST 48 (2015) 21–83.
  * Ju et al. [2016] Y. Ju, J. K. Lefkowitz, C. B. Reuter, S. H. Won, X. Yang, S. Yang, W. Sun, Z. Jiang, Q. Chen,  Plasma assisted low temperature combustion,  PLASMA CHEM PLASMA P 36 (2016) 85–105.
  * Starikovskiy and Aleksandrov [2013] A. Starikovskiy, N. Aleksandrov,  Plasma-assisted ignition and combustion,  PROG ENERG COMBUST 39 (2013) 61–110.
  * Mao and Chen [2018] X. Mao, Q. Chen,  Effects of vibrational excitation on nanosecond discharge enhanced methane–air ignition,  AIAA J 56 (2018) 4312–4320.
  * Mao et al. [2019a] X. Mao, Q. Chen, C. Guo,  Methane pyrolysis with n2/ar/he diluents in a repetitively-pulsed nanosecond discharge: Kinetics development for plasma assisted combustion and fuel reforming,  ENERG CONVERS MANAGE 200 (2019a) 112018.
  * Mao et al. [2019b] X. Mao, A. Rousso, Q. Chen, Y. Ju,  Numerical modeling of ignition enhancement of ch4/o2/he mixtures using a hybrid repetitive nanosecond and dc discharge,  P COMBUST INST 37 (2019b) 5545–5552.
  * Lefkowitz et al. [2015] J. K. Lefkowitz, P. Guo, A. Rousso, Y. Ju,  Species and temperature measurements of methane oxidation in a nanosecond repetitively pulsed discharge,  PHILOS T ROY SOC A 373 (2015) 20140333.
  * Sun et al. [2022] J. Sun, Q. Chen, B. Zhao, C. Guo, J. Liu, M. Zhang, D. Li,  Temperature-dependent ion chemistry in nanosecond discharge plasma-assisted ch4 oxidation,  J PHYS D APPL PHYS 55 (2022) 135203.
  * Welz et al. [2015] O. Welz, M. P. Burke, I. O. Antonov, C. F. Goldsmith, J. D. Savee, D. L. Osborn, C. A. Taatjes, S. J. Klippenstein, L. Sheps,  New insights into low-temperature oxidation of propane from synchrotron photoionization mass spectrometry and multiscale informatics modeling,  J PHYS CHEM A 119 (2015) 7116–7129.
  * Merchant et al. [2015] S. S. Merchant, C. F. Goldsmith, A. G. Vandeputte, M. P. Burke, S. J. Klippenstein, W. H. Green,  Understanding low-temperature first-stage ignition delay: Propane,  COMBUST FLAME 162 (2015) 3658–3673.
  * Adamovich et al. [2015] I. V. Adamovich, T. Li, W. R. Lempert,  Kinetic mechanism of molecular energy transfer and chemical reactions in low-temperature air-fuel plasmas,  PHILOS T ROY SOC A 373 (2015) 20140336.
  * Eckert et al. [2018] Z. Eckert, N. Tsolas, K. Togai, A. Chernukho, R. A. Yetter, I. V. Adamovich,  Kinetics of plasma-assisted oxidation of highly diluted hydrocarbon mixtures excited by a repetitive nanosecond pulse discharge,  J PHYS D APPL PHYS 51 (2018) 374002.
  * Chen et al. [2023] H. Chen, R. Zhang, H. Liao, F. Zhang, J. Yang, B. Yang,  Kinetic insights into plasma-assisted low-temperature oxidation of propane with synchrotron photoionization mass spectrometry,  P COMBUST INST 39 (2023) 5499–5509.
  * Ban et al. [2023] Y. Ban, S. Zhong, J. Zhu, F. Zhang,  Effects of non-equilibrium plasma and equilibrium discharge on low-temperature combustion in lean propane/air mixtures,  Fuel 339 (2023) 127353\.
  * Wu et al. [2011] L. Wu, J. Lane, N. Cernansky, D. Miller, A. Fridman, A. Y. Starikovskiy,  Plasma-assisted ignition below self-ignition threshold in methane, ethane, propane and butane-air mixtures,  P COMBUST INST 33 (2011) 3219–3224.
  * Akishev et al. [2022] Y. Akishev, T. Alekseeva, V. Karalnik, A. Petryakov,  On the slow ionization waves forming the breakdown in a long capillary tube with helium at low pressure,  J PHYS D APPL PHYS 55 (2022) 145202.
  * Klochko et al. [2015] A. Klochko, A. Salmon, J. Lemainque, N. Popov, S. Starikovskaia,  Capillary nanosecond discharges as a tool for the measurement of quenching coefficients at high specific energy deposition,  in: 53rd AIAA Aerospace Sciences Meeting, 2015, p. 0937.
  * Chng et al. [2019] T. L. Chng, I. Orel, S. Starikovskaia, I. Adamovich,  Electric field induced second harmonic (e-fish) generation for characterization of fast ionization wave discharges at moderate and low pressures,  PLASMA SOURCES SCI T 28 (2019) 045004.
  * Klochko et al. [2014] A. Klochko, A. Salmon, N. Popov, J.-P. Booth, M. J. Kushner, Z. Xiong, S. Starikovskaya,  Experimental and numerical study of fast gas heating and o atom production in a capillary nanosecond discharge,  in: 52nd Aerospace Sciences Meeting, 2014, p. 1030.
  * Lepikhin et al. [2016] N. Lepikhin, A. Klochko, N. Popov, S. Starikovskaia,  Long-lived plasma and fast quenching of n2 (c3$\pi$u) by electrons in the afterglow of a nanosecond capillary discharge in nitrogen,  PLASMA SOURCES SCI T 25 (2016) 045003.
  * Chen et al. [2021] X. Chen, Y. Zhu, Y. Wu, J. Hao, X. Ma, P. Lu,  Modeling of fast ionization waves in pure nitrogen at moderate pressure,  PLASMA SOURCES SCI T 30 (2021) 065002.
  * Zhu et al. [2020] Y. Zhu, S. M. Starikovskaia, N. Y. Babaeva, M. J. Kushner,  Scaling of pulsed nanosecond capillary plasmas at different specific energy deposition,  PLASMA SOURCES SCI T 29 (2020) 125006.
  * Pancheshnyi et al. [2008] S. Pancheshnyi, B. Eismann, G. Hagelaar, L. Pitchford,  Computer code zdplaskin,  University of Toulouse, LAPLACE, CNRS-UPS-INP, Toulouse, France (2008). URL: http://www.zdplaskin.laplace.univ-tlse.fr.
  * Luque [2020] A. Luque,  Qtplaskin,  the Instituto de Astrofísica de Andalucía (IAA) (2020). URL: https://github.com/aluque/qtplaskin/.
  * Markosyan et al. [2014] A. H. Markosyan, A. Luque, F. J. Gordillo-Vázquez, U. Ebert,  Pumpkin: A tool to find principal pathways in plasma chemical models,  Computer Physics Communications 185 (2014) 2697–2702.
  * Qiu et al. [2023] Y. Qiu, Y. Zhu, Y. Wu, N. Zhao, Z. Li, M. Hao, B. Zhang, D. Pan,  Numerical investigation of the hybrid pulse–dc plasma assisted ignition and nox emission of nh3/n2/o2 mixture,  COMBUST FLAME 258 (2023) 113078.
  * Hagelaar and Pitchford [2005] G. Hagelaar, L. C. Pitchford,  Solving the boltzmann equation to obtain electron transport coefficients and rate coefficients for fluid models,  PLASMA SOURCES SCI T 14 (2005) 722.
  * Hindmarsh [1983] A. C. Hindmarsh,  Odepack, a systemized collection of ode solvers,  Scientific computing (1983).
  * Guglielmi and Hairer [2001] N. Guglielmi, E. Hairer,  Implementing radau iia methods for stiff delay differential equations,  Computing 67 (2001) 1–12.
  * Yang and Tao [2006] S. Yang, W. Tao, Chuan re xue[Heat transfer theory], Higher Education Press, 2006\.
  * Reuter et al. [2018] C. B. Reuter, R. Zhang, O. R. Yehia, Y. Rezgui, Y. Ju,  Counterflow flame experiments and chemical kinetic modeling of dimethyl ether/methane mixtures,  COMBUST FLAME 196 (2018) 1–10.
  * Zhou et al. [2018] C.-W. Zhou, Y. Li, U. Burke, C. Banyon, K. P. Somers, S. Ding, S. Khan, J. W. Hargis, T. Sikes, O. Mathieu, et al.,  An experimental and chemical kinetic modeling study of 1, 3-butadiene combustion: Ignition delay time and laminar flame speed measurements,  COMBUST FLAME 197 (2018) 423–438.
  * Ramalingam et al. [2021] A. Ramalingam, S. Panigrahy, Y. Fenard, H. Curran, K. A. Heufer,  A chemical kinetic perspective on the low-temperature oxidation of propane/propene mixtures through experiments and kinetic analyses,  COMBUST FLAME 223 (2021) 361–375.
  * Janev and Reiter [2002] R. K. Janev, D. Reiter, Collision processes of hydrocarbon species in hydrogen plasmas: II. The ethane and propane families, Forschungszentrum, Zentralbibliothek, 2002\.
  * Zhu et al. [2018] Y. Zhu, N. D. Lepikhin, I. S. Orel, A. Salmon, A. V. Klochko, S. M. Starikovskaia,  Optical actinometry of o-atoms in pulsed nanosecond capillary discharge: peculiarities of kinetics at high specific deposited energy,  PLASMA SOURCES SCI T 27 (2018) 075020.
  * Klippenstein et al. [2011] S. J. Klippenstein, L. B. Harding, M. J. Davis, A. S. Tomlin, R. T. Skodje,  Uncertainty driven theoretical kinetics studies for ch3oh ignition: Ho2+ ch3oh and o2+ ch3oh,  P COMBUST INST 33 (2011) 351–357.
  * Bedjanian and Morin [2017] Y. Bedjanian, J. Morin,  Reaction of o (3p) with c3h6: Yield of the reaction products as a function of temperature,  J PHYS CHEM A 121 (2017) 1553–1562.
  * Leonori et al. [2015] F. Leonori, N. Balucani, V. Nevrly, A. Bergeat, S. Falcinelli, G. Vanuzzo, P. Casavecchia, C. Cavallotti,  Experimental and theoretical studies on the dynamics of the o (3p)+ propene reaction: primary products, branching ratios, and role of intersystem crossing,  J PHYS CHEM C 119 (2015) 14632–14652.
  * Tsang [1987] W. Tsang,  Chemical kinetic data base for combustion chemistry. part 2. methanol,  J PHYS CHEM REF DATA 16 (1987) 471–508.
  * Atkinson et al. [2006] R. Atkinson, D. Baulch, R. Cox, J. Crowley, R. Hampson, R. Hynes, M. Jenkin, M. Rossi, J. Troe, I. Subcommittee,  Evaluated kinetic and photochemical data for atmospheric chemistry: Volume ii–gas phase reactions of organic species,  ATMOS CHEM PHYS 6 (2006) 3625–4055.


# Constructions for Nonadaptive
Tropical Group Testing

Nicholas Kwan and Lele Wang  University of British Columbia, Vancouver, BC
V6T1Z4, Canada
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

PCR testing is an invaluable diagnostic tool that has most recently seen
widespread use during the COVID-19 pandemic. A recent work by Wang, Gabrys and
Vardy proposed tropical codes as a model for group PCR testing. For a known
but arbitrary number of infected persons, a sufficient condition on the
underlying block design of a zero-error tropical code, called double
disjunction, is proposed. Despite this, the parameters for which the
construction of doubly disjunct block designs is known to exist are very
limited. In this paper, we define probabilistic tropical codes and consider
random block designs that are doubly disjunct with high probability. We also
provide a deterministic construction for a doubly disjunct block design given
a disjunct block design. We show that for certain choices of parameters, our
probabilistic construction has vanishing error. Our constructions, combined
with existing methods, give us three different ways to construct tropical
codes. We compare the number of tests required by each, and bounds on the
error.

##  I Introduction

Polymerase Chain Reaction (PCR) is a method by which genetic material may be
rapidly amplified. PCR has a wide range of applications, perhaps most recently
in reliable testing for the virus causing COVID-19. PCR testing is typically
done through thermal cycling, where the viral load of a test containing some
specimens is repeatedly doubled in concentration. Each doubling event is known
as a cycle, and the load is doubled until it reaches a detectable threshold or
fails to reach this threshold within a preset number of cycles. The cycle
number at which this threshold is reached is called the Ct value. The Ct value
of a test is thus

$\text{Ct value}=\lfloor-\log_{2}(\text{viral load})+\text{constant}\rfloor,$

and is a measure of the degree of infection of a test, with smaller values
indicating higher degrees of infection. If a test fails to reach the threshold
after a preset number of cycles, then it is not infected and its Ct-value is
$\infty$.

Consider some finite set $\mathcal{S}$ of people with an infected subset
$\mathcal{D}\subset\mathcal{S}$. Group testing deals with the problem of
identifying $\mathcal{D}$ by performing tests on subsets of $\mathcal{S}$
whose size may be greater than one. By nonadaptive group testing, we consider
a setting where all group tests are given at once. Group testing can be
divided into combinatorial and probabilistic group testing. Combinatorial
group testing identifies $\mathcal{D}$ with zero error, while probabilistic
group testing identifies $\mathcal{D}$ with high probability. Traditionally,
group testing takes place in the binary setting so that a test on some subset
$\mathcal{U}\subset\mathcal{S}$ is positive when $\mathcal{U}\cap\mathcal{D}$
is nonempty and negative otherwise. On the other hand, a PCR test on a subset
$\mathcal{U}$ returns some aggregate Ct-value, which tells us more than
whether or not $\mathcal{U}\cap\mathcal{D}$ is nonempty. For instance, when
two specimens with Ct values $x$ and $y$ are placed into a PCR test at the
same time, the resulting Ct value is $z$ satisfying

$2^{-x}+2^{-y}=2^{-z}.$

Due to exponential decay, we have $z\approx\min\\{x,y\\}$. In performing a
group PCR test, we are also able to add a specimen to a test after a certain
number of cycles have elapsed. The operation of adding a specimen to a test
after some $\delta$ cycles is known as delaying, and we say that the specimen
has been delayed by $\delta$ cycles. When a specimen with Ct value $x$ is
delayed by $\delta$ cycles, the number of cycles it takes to reach the same
detectable threshold becomes $x+\delta$.

Based on the above properties of PCR testing, Wang, Gabrys, and Vardy proposed
the tropical semiring as a model for group PCR testing [1]. The tropical
semiring $(\mathbb{R}\cup\\{\infty\\},\oplus,\odot)$ consists of the real
numbers and infinity equipped with the following operations:

$\displaystyle x\oplus y$ $\displaystyle:=\min(x,y),$ $\displaystyle x\odot y$
$\displaystyle:=x+y.$

The additive identity is $\infty$ and the multiplicative identity is $0$. In
this model, for two specimens $x$ and $y$ (we have identified the specimens
with their individual Ct-values), the Ct-value
$(x\odot\delta_{x})\oplus(y\odot\delta_{y})=\min(x+\delta_{x},y+\delta_{y})$
is given by a test containing both $x$ and $y$, delayed each by $\delta_{x}$
and $\delta_{y}$ respectively. In this way, $\odot$ indicates delay and
$\oplus$ indicates presence. This model works because viral loads vary widely,
and the Ct-value of a test with two specimens with differing viral loads is
dominated by the specimen with the higher viral load due to the exponential
nature of the test. Since Ct values and delays only take nonnegative integer
values, all the analysis done in this paper is over the subsemiring
$(\mathbb{N}\cup\\{0,\infty\\},\oplus,\odot)$.

In general, a test involving specimens $x_{1},\ldots,x_{k}$ with delays
$\delta_{1},\ldots,\delta_{k}$ will return the Ct-value
$\bigoplus_{i=1}^{k}(x_{i}\odot\delta_{i})$. In this way, we can define tests
using tropical matrix multiplication. Consider a set of $n$ individuals, and
let
$x=[x_{1},x_{2},\ldots,x_{n}]^{\intercal}\in(\mathbb{N}\cup\\{0,\infty\\})^{n}$
be a list of their Ct values. A set of $t$ tests on these individuals can be
formally written as a matrix of delay values
$S\in(\mathbb{N}\cup\\{0,\infty\\})^{t\times n}$, known as the schedule
matrix. If $S_{ij}=\infty$, then the $i^{\text{th}}$ test did not involve the
$j^{\text{th}}$ individual. The results of the tests are then given by $S\odot
x$, where here, $\odot$ denotes tropical matrix multiplication; that is,

$S\odot x=\begin{bmatrix}\bigoplus_{j=1}^{N}(S_{1j}\odot x_{j})\\\ \vdots\\\
\bigoplus_{j=1}^{N}(S_{tj}\odot x_{j})\end{bmatrix}.$

If, for any distinct $x,y\in(\mathbb{N}\cup\\{0,\infty\\})^{n}$ with at most
$d$ finite values, $S\odot x\neq S\odot y$, then our schedule matrix is a
valid $(t,n,d)$-tropical code in the combinatorial setting. The code is said
to have maximum delay $l$, where $l$ is the largest finite element in $S$. In
a real-world setting, $l$ is typically held at 40 cycles. For general $d$,
zero-error tropical codes are hard to construct.

Existing work in nonadaptive group testing considers the construction of
binary matrices satisfying certain properties, such as disjunction and
separability, that correspond to binary group testing schemes. These
properties are necessary, but not sufficient, for binary matrices that
correspond to the schedule matrix of a tropical code. Furthermore, the
sufficient doubly disjunct property given in definition 1 by [1] is new, and
to our knowledge, existing constructions do not exist for doubly disjunct
matrices. In this paper, we propose the construction of _random_ schedule
matrices that are doubly disjunct with high probability, as well as the
construction of doubly disjunct matrices from disjunct matrices.

Notation. Throughout we use calligraphic letters such as $\mathcal{K}$ to
refer to sets and block designs. We use boldfaced letters for random matrices
$\bm{S}$ and random vectors $\bm{s}$ when they are capitalized and in
lowercase, respectively. $\operatorname{\textsf{P}}(A)$ refers to the
probability of the event $A$, and we define
$p_{X}(x)\triangleq\operatorname{\textsf{P}}(X=x)$. Plain capital letters are
used for matrices and blocks, and $\vec{p}$ is used for vectors.

##  II Summary of Existing Results

This paper relies heavily on the results given by [1], which was the first
paper to consider using the tropical semiring as a model for group PCR
testing, and produced results on tropical codes in the combinatorial setting.

For some $t\in\mathbb{N}$, a block design on $t$ vertices is a subset
$\mathcal{F}\subset 2^{[t]}$ of the power set $2^{[t]}$ satisfying certain
properties with respect to unions, intersections, and differences. Each
element $B\in\mathcal{F}$ is called a block. The set operation properties
satisfied by the blocks vary between block designs. An example of a block
design on $7$ vertices is the Fano plane
$\mathcal{F}=\\{\\{1,2,4\\},\\{1,3,7\\},\\{1,5,6\\},\\{2,3,5\\},\\{2,6,7\\},\\{3,4,6\\},$
$\\{4,5,7\\}\\}$. Note that the intersection of any two blocks has exactly one
element, and that each block has size three. Furthermore, each element of
$[7]$ is contained in exactly three blocks.

We may associate a block design $\mathcal{F}$ on $T$ vertices with
$|\mathcal{F}|=N$ with an incidence matrix $M$, a $T\times N$ binary matrix
where each column corresponds to a block. Then

$M_{ij}=\begin{cases}0&\text{if $i\notin B_{j}$ for the block
$B_{j}\in\mathcal{F},$}\\\ 1&\text{if $i\in B_{j}$ for the block
$B_{j}\in\mathcal{F}.$}\end{cases}$

Note that this representation is unique up to permutation of the columns $M$,
whereas the permutation of the rows of $M$ results in an equivalent block
design. For instance, the incidence matrix of the previously described Fano
plane is

$\begin{bmatrix}1&0&0&0&1&0&1\\\ 1&1&0&0&0&1&0\\\ 0&1&1&0&0&0&1\\\
1&0&1&1&0&0&0\\\ 0&1&0&1&1&0&0\\\ 0&0&1&0&1&1&0\\\
0&0&0&1&0&1&1\end{bmatrix}.$

We may use the incidence matrices of suitable block designs to construct
tropical codes by replacing each $0$ with $\infty$ and each $1$ with a
suitable finite delay value. Conversely, every tropical code has an associated
block design.

######  Definition 1 (Doubly Disjunct Block Design).

A block design $\mathcal{F}$ with $n$ blocks on $t$ vertices is said to be
$d$-doubly disjunct if, for any distinct blocks
$Z,B_{1},\ldots,B_{d}\in\mathcal{F}$,
$\left|Z\setminus\left(\bigcup_{i=1}^{d}B_{i}\right)\right|\geq 2.$

If instead $\left|Z\setminus\left(\bigcup_{i=1}^{d}B_{i}\right)\right|\geq 1$,
then we call $\mathcal{F}$ $d$-disjunct.

######  Remark.

Clearly if $\mathcal{F}$ is a $d$-doubly disjunct block design, it is also
$d$-disjunct.

The Fano plane is a 1-doubly disjunct block design on $7$ vertices. A
construction for $1$-doubly disjunct block design on some arbitrary number of
vertices $t$ is given by Graham and Sloane [2]. Let $\mathcal{S}^{t}_{w}$ be
the set of all binary vectors of length $t$ and Hamming weight $w$, consider
the map

$\displaystyle f:\mathcal{S}_{w}^{t}$
$\displaystyle\longrightarrow\mathbb{Z}/t\mathbb{Z},$ $\displaystyle\vec{v}$
$\displaystyle\longmapsto\left(\sum_{i=1}^{n}iv_{i}\right)\mod t,$

taking a binary vector $v$ to the sum of the indices of its ones, modulo $t$.
Then the preimage of any residue class will be a $1$-doubly disjunct block
design. At least one of these preimages will have cardinality larger than
$\frac{1}{t}\binom{t}{w}$.

######  Theorem 1 (Theorem 31 in [1]).

Let $\mathcal{F}$ be a $(d-1)$-doubly-disjunct block design with $n$ blocks on
$t$ vertices. Let $M$ be the incidence matrix of $\mathcal{F}$. Then there
exists a $(t,n,d)$-tropical code whose associated block design is
$\mathcal{F}$. That is, the finite entries of the schedule matrix $S$ coincide
exactly with the nonzero entries of $M$.

In the case where $d=2$ and each block $B\in\mathcal{F}$ has size at least
three, there exists a decoding algorithm, given by Theorem 17 in [1], whose
maximum delay is a prime bounded above by $n$. Otherwise, for $d>2$, the
decoding algorithm given in Theorem 31 in [1] uses a maximum delay of
$2^{n(t+1)}$, with the schedule matrix

$S_{ij}=\begin{cases}2^{i+jt},&\text{ if }M_{ij}=1,\\\ \infty,&\text{ if
}M_{ij}=0.\end{cases}$

although the authors conjecture that there is a polynomial delay code using
the same underlying block design. The given decoding algorithm considers sets
of potentially infected persons of size $d$. We may construct a bipartite
graph using a plausible set of $d$ infected persons and the actual set of $d$
infected persons as our two sets of vertices. The double disjunctness of
$\mathcal{F}$ guarantees that this graph contains an even cycle. The
alternating sum along this cycle will correspond to the sum of delays for the
actual set of infected persons.

To the best of our knowledge, there is no construction for a $(d-1)$-doubly
disjunct block design for $d>2$ with arbitrary parameters. Existing
constructions for $d$-disjunct matrices, which are suitable for binary group
testing, exist. Combinatorial constructions exist with $t=O(d^{2}\log n)$ [3],
while Monte-Carlo constructions, like our own, exist for
$O\left(d^{2}\min\\{\log n,(\log_{t}n)^{2}\\}\right)$ [4]. Probabilistic
binary group testing schemes exist with $t=\Theta(d\log n)$ [5].

##  III Result I: A Deterministic Construction from Disjunct Block Designs

Recall from Definition 1 what it means for a block design $\mathcal{F}$ to be
disjunct and doubly disjunct. Our first result provides a simple way to
construct $d$-doubly disjunct block designs from $d$-disjunct block designs by
doubling the number of vertices. These can then be used to produce
$(t,n,d)$-tropical codes by Theorem 1. Methods for constructing $d$-disjunct
block designs are given in [5, 4, 3].

######  Theorem 2.

Let $\mathcal{F}\subset 2^{[t]}$ be a $d$-disjunct block design on $t$
vertices. Then the image $\mathcal{G}$ of the function

$\displaystyle f:\mathcal{F}$ $\displaystyle\longrightarrow 2^{[2t]}$
$\displaystyle B$ $\displaystyle\longmapsto\\{2x:x\in B\\}\cup\\{2x-1:x\in
B\\},$

is a $d$-doubly disjunct block design on $2t$ vertices.

###### Proof:

Let $B_{0},\ldots,B_{d}$ be any $d+1$ distinct blocks in $\mathcal{G}$, and
for each $i\in\\{0,\ldots,d\\}$ let $A_{i}$ be the unique element of
$\mathcal{F}$ mapping to $B_{i}$. Since $\mathcal{F}$ is disjunct,
$\left|A_{0}\setminus\left(\bigcup_{i=1}^{d}A_{i}\right)\right|\geq 1$. For
each block, the maps $x\mapsto 2x$ and $x\mapsto 2x-1$ are injective, and the
images of the two maps are disjoint. Let $a\in
A_{0}\setminus\left(\bigcup_{i=1}^{d}A_{i}\right)$. Then note that $2a\in
B_{0}$ and $2a-1\in B_{0}$. If $2a\in\bigcup_{i=1}^{d}B_{i}$, then $2a\in
B_{i}$ for some $i\geq 1$, but this would imply that $a\in A_{i}$ for some
$i\geq 1$, which cannot be the case. Similarly, if
$2a-1\in\bigcup_{i=1}^{d}B_{i}$, then $2a-1\in B_{i}$ for some $i\geq 1$, but
this would imply that $a\in A_{i}$ for some $i\geq 1$, which cannot be the
case. Hence, $2a\in B_{0}\setminus\left(\bigcup_{i=1}^{d}B_{i}\right)$ and
$2a-1\in B_{0}\setminus\left(\bigcup_{i=1}^{d}B_{i}\right)$, so
$\left|B_{0}\setminus\left(\bigcup_{i=1}^{d}B_{i}\right)\right|\geq 2$ and
$\mathcal{G}$ is $d$-doubly disjunct. ∎

######  Remark.

This construction extends to the construction of $n$-fold disjunct block
designs, but it is not optimal. For instance, consider the Fano plane as
described in section II. The Fano plane is an “optimal” 2-disjunctive block
design on 7 vertices in the sense that no blocks can be added to it while
preserving 2-disjunctiveness, but applying Theorem 2 to the Fano plane and
adding the block $\\{2,4,6,8,10,12,14\\}$ still yields a 2-doubly disjunct
block design on 14 vertices.

##  IV Result II: A Probabilistic Construction

In this section we define tropical codes in a probabilistic setting, and
construct random matrices that are doubly disjunct with high probability.

###  IV-A Definition

Consider a _random_ $t\times n$ schedule matrix $\mathbf{S}$. We define the
probability of error as

$\displaystyle
P_{e}^{n,d}(\bm{S})\triangleq\operatorname{\textsf{P}}\left({\bm{S}}\odot\vec{x}={\bm{S}}\odot\vec{y}\text{
for some distinct }\vec{x},\vec{y}\in\mathcal{R}_{n,d}\right),$

where
$\mathcal{R}_{n,d}\triangleq\\{\vec{x}\in(\mathbb{N}\cup\\{0,\infty\\})^{n}:x\text{
contains at most $d$}\\\ \text{finite values}\\}$. In other words,
$P_{e}^{n,d}(\bm{S})$ is the probability that a group testing scheme with
schedule matrix $\bm{S}$ confuses $x$ and $y$. To formalize the definition of
a tropical code in a probabilistic setting, we introduce a parameter
$\epsilon$ to bound the probability of error $P_{e}^{n,d}(\bm{S})$.

######  Definition 2.

A random $t\times n$ schedule matrix $\bm{S}$ is a $(t,n,d,\epsilon)$-tropical
code if $P_{e}^{n,d}\leq\epsilon$.

A $(t,n,d,0)$-tropical code then corresponds to the definition of a
$(t,n,d)$-tropical code given in [1]. This definition naturally motivates us
to find relationships between the parameters $t,n,d$ and $\epsilon$.

###  IV-B Construction

In this section we find parameters $(t,n,d)$ for which we are able to
construct doubly disjunct matrices with high probability, then compute
$\epsilon$ for our construction, based on its probability of generating a
matrix that is not doubly disjunct.

We show the existence of a joint distribution $p_{\bm{s}}(\vec{v})$ for the
first row of delays $\bm{s}$ in $\bm{S}$. We then generate the remaining rows
of $\bm{S}$ i.i.d. according to the same distribution $p_{\bm{s}}(\vec{v})$.

We produce sufficient properties on $p_{\bm{s}}(\vec{v})$ for $\bm{S}$ to be
doubly disjunct with high probability. Since $p_{\bm{s}}(\vec{v})$ is a
probability distribution, we require

$\sum_{\vec{v}\in\\{0,1\\}^{n}}p_{\bm{s}}(\vec{v})=1,$ (1)

and

$p_{\bm{s}}(\vec{v})\geq 0\text{ for all }\vec{v}\in\\{0,1\\}^{n}.$ (2)

For $\bm{S}$ to be $d-1$-doubly disjunct on average, we require, for every
distinct $i_{1},i_{2},\ldots,i_{d}\in[n]$

$\displaystyle\operatorname{\textsf{P}}\begin{pmatrix}\bm{s}_{i_{1}}=1\text{
and}\\\
\bm{s}_{i_{2}}=\cdots=\bm{s}_{i_{d}}=0\end{pmatrix}=\sum_{\begin{subarray}{c}\vec{v}\in\\{0,1\\}^{n}\\\
v_{i_{1}}=1,\\\
v_{i_{2}}=\cdots=v_{i_{d}}=0\end{subarray}}p_{\bm{s}}(\vec{v})=\frac{2+\Delta}{t}.$
(3)

Additionally, to ensure that no test is empty, we set $P(\bm{s}=\vec{0})=0$.
Conditions (1), (2), and (3) reduce our problem to proving the existence of
non-negative solutions to a linear system. To simplify our analysis, we
restrict to the case where rows whose Hamming weights are equal have equal
probability of being generated. For instance, if $n=3$, then we let
$p_{\bm{s}}(100)=p_{\bm{s}}(010)=p_{\bm{s}}(001)$ and
$p_{\bm{s}}(110)=p_{\bm{s}}(101)=p_{\bm{s}}(011)$. This allows us to write our
first conditions (1) and (3) as linear combinations of the probabilities that
$\bm{s}$ has hamming weight $w\leq n$, rather than the probabilities of
$p_{\bm{s}}(\vec{v})$ for each $\vec{v}\in\\{0,1\\}^{n}$. Let $\vec{p}$ be the
vector whose $w^{\text{th}}$ element corresponds to the probability that
$\vec{s}$ has hamming weight $w$. Conditions (1) and (3) can then be written
as the following linear system, for which $\vec{p}$ is a solution

$\begin{bmatrix}\binom{n}{1}&\binom{n}{2}&\cdots&\binom{n}{n-d+1}&\;\;\;\cdots&\binom{n}{n}\\\\[5.0pt]
\binom{n-d}{0}&\binom{n-d}{1}&\cdots&\binom{n-d}{n-d}&0\;\;\cdots&0\end{bmatrix}\vec{p}=\begin{bmatrix}1\\\
\frac{2+\Delta}{t}\end{bmatrix}.$ (4)

The first row of the matrix in equation (4) corresponds to the condition given
by (1); its $w^{\text{th}}$ component is the number of length $n$ binary
vectors with weight $w$. The second row corresponds to the condition given by
(3); its $w^{\text{th}}$ component is the number of weight $w$ length $n$
binary vectors $\bm{s}$ with $\bm{s}_{i_{1}}=1$ and
$\bm{s}_{i_{2}}=\ldots=\bm{s}_{i_{d}}=0$ for any distinct
$i_{1},\ldots,i_{d}\in[n]$.

The following lemma provides conditions on the parameters such that an
elementwise nonnegative solution to the linear system (4) exists.

######  Lemma 1.

For a general pool of $n$ patients with a maximum of $d$ infected, and any
$\Delta>0$, if

$\frac{\binom{n-d}{\lfloor\frac{n}{d}\rfloor-1}}{\binom{n}{\lfloor\frac{n}{d}\rfloor}}>\frac{2+\Delta}{t},$

then equation (4) has an elementwise nonnegative solution.

We will use Farkas’ Lemma to establish Lemma 1.

######  Lemma 2 (Farkas’ Lemma, [6, p. 263]).

Let $A\in\mathbb{R}^{m\times n}$ and $\vec{b}\in\mathbb{R}^{m}$. Then exactly
one of the following is true:

  1. 1.

There exists $\vec{p}\in\mathbb{R}^{n}$ such that $A\vec{p}=\vec{b}$ and
$\vec{p}\geq 0$ elementwise

  2. 2.

There exists $\vec{y}\in\mathbb{R}^{m}$ such that $A^{\intercal}\vec{y}\geq 0$
elementwise and $\vec{b}^{\intercal}\vec{y}<0$

###### Proof:

If no elementwise nonnegative vector $\vec{p}$ exists satisfying (4), then by
the second condition of Farkas’ lemma, there exists some
$\vec{y}=[y_{1}\;y_{2}]^{\intercal}$ such that $y_{1}\geq 0$,
$y_{1}+\frac{2+\Delta}{t}y_{2}<0$, and
$\binom{n}{w}y_{1}+\binom{n-d}{w-1}y_{2}\geq 0$ for all $w=1,2,\ldots,n-d+1$.
Together, these imply

$\displaystyle 0$ $\displaystyle>y_{1}+\frac{2+\Delta}{t}y_{2}$
$\displaystyle\geq\frac{-\binom{n-d}{w-1}y_{2}}{\binom{n}{w}}+\frac{2+\Delta}{t}y_{2}$
$\displaystyle=y_{2}\left(\frac{2+\Delta}{t}-\frac{\binom{n-d}{w-1}}{\binom{n}{w}}\right),$

for all $w=1,2,\ldots,n-d+1$. Note that, since $y_{1}\geq 0$ and
$y_{1}+\frac{2+\Delta}{t}y_{2}<0$, we must have $y_{2}<0$. This means

$\frac{\binom{n-d}{w-1}}{\binom{n}{w}}<\frac{2+\Delta}{t}.$

We compute the maximum of
$f(w)\triangleq\frac{\binom{n-d}{w-1}}{\binom{n}{w}}$ over $w\in[d]$. Note
that

$\displaystyle f(w)$
$\displaystyle=\frac{w\prod_{i=0}^{w-2}(n-d-i)}{\prod_{i=0}^{w-1}(n-i)},$

so that

$\displaystyle f(w+1)$
$\displaystyle=\frac{(w+1)(n-d-w+1)\prod_{i=0}^{w-2}(n-d-i)}{(n-w)\prod_{i=0}^{w-1}(n-i)}.$

We examine the terms outside of the products in the numerator and denominator,
which differ between $f(w)$ and $f(w+1)$. For $f(w+1)$, we have

$\displaystyle\frac{(w+1)(n-d-w+1)}{n-w}$
$\displaystyle=w\left(1+\frac{1}{w}\right)\left(1-\frac{d-1}{n-w}\right).$

If $w\geq\lfloor n/d\rfloor$ is an integer, then $wd\geq n$ so that

$\displaystyle w(d-1)$ $\displaystyle\geq n-w$ $\displaystyle\frac{d-1}{n-w}$
$\displaystyle\geq\frac{1}{w},$

and $f(w+1)\leq f(w)$. If $w<\lfloor n/d\rfloor$ is an integer, then
$w\leq\frac{n}{d}-1$ so that

$\displaystyle n-w$ $\displaystyle\geq d(w+1)-w$ $\displaystyle>(d-1)(w+1)$
$\displaystyle\frac{1}{w}$
$\displaystyle>\left(\frac{d-1}{n-w}\right)\left(1+\frac{1}{w}\right),$

and $f(w+1)>f(w)$. As such, $f$ is maximized when $w=\lfloor n/d\rfloor$.
Hence, if

$\frac{\binom{n-d}{\lfloor\frac{n}{d}\rfloor-1}}{\binom{n}{\lfloor\frac{n}{d}\rfloor}}>\frac{2+\Delta}{t},$

then there exists some elementwise non-negative vector $\vec{p}$ satisfying
(4). Thus, there exists a valid probability distribution on $\mathbf{s}$,
using which we can generate a $(t,n,d)$-tropical code with high probability. ∎

###  IV-C Performance of the probabilistic construction

In this section, we establish conditions for which the proposed construction
provides us with a valid tropical code.

######  Theorem 3.

The construction given in Section IV-B is a
$\left(t,n,d,\epsilon\right)$-tropical code if

$t>2e\left(d^{2}\log\left(\frac{ne}{d}\right)+d\log\left(\frac{d}{\epsilon}\right)+2d\right).$

We will use Chernoff bound to establish Theorem 3.

######  Lemma 3 (Chernoff bound [7]).

Let $X=\sum_{k=1}^{t}X_{k}$ where all $X_{k}\sim\operatorname{Bern}(p)$ are
independent. Let $\mu=\mathbb{E}(X)=tp$. Then

$\operatorname{\textsf{P}}(X\leq(1-\delta)\mu)\leq e^{-\mu\delta^{2}/2},$

for all $0<\delta<1.$

###### Proof:

For a given $\epsilon\in(0,1)$, choose
$\Delta=\left(2+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)\right)$. For
some $d$ fixed distinct blocks indexed $i_{1},\ldots,i_{d}\in[n]$ let
$X_{k}=\mathbf{1}_{\\{{\bm{S}}_{k,i_{1}}=1,{\bm{s}}_{k,i_{2}}=0,\ldots{\bm{s}}_{k,i_{d}}=0\\}}\sim\operatorname{Bern}(\frac{4+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)}{t})$
be the random variable indicating whether or not the $i_{1}^{\text{th}}$ entry
of the $k^{\text{th}}$ row is $1$ and the $i_{2},\ldots,i_{d}$ entries of the
$k^{\text{th}}$ row are 0. Then the probability of error for the chosen subset
of blocks is $\operatorname{\textsf{P}}\left(\sum_{k\in[t]}X_{k}<2\right)$,
since the sum $\sum_{k\in[t]}X_{k}$ gives the size of the set
$\left|B_{i_{1}}\setminus\left(\bigcup_{j=2}^{d}B_{i_{j}}\right)\right|$. Then
note that

$\displaystyle\operatorname{\textsf{P}}\Big{(}\sum_{k\in[t]}X_{k}<2\Big{)}$
$\displaystyle=\operatorname{\textsf{P}}\left(\sum_{k\in[t]}X_{k}<\Big{(}4+2\log\Big{(}\tfrac{d}{\epsilon}\tbinom{n}{d}\Big{)}\Big{)}\left(1-\tfrac{2+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)}{4+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)}\right)\right).$

We may then use the Chernoff bound (lemma 3) with
$\mu=4+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)$ and
$\delta=\frac{2+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)}{4+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)}$
to get:

$\displaystyle\operatorname{\textsf{P}}\Big{(}\sum_{k\in[t]}X_{k}<2\Big{)}$
$\displaystyle\leq\exp\left(-\left(4+2\log\left(\tfrac{d}{\epsilon}\tbinom{n}{d}\right)\right)\tfrac{\left(2+2\log\left(\tfrac{d}{\epsilon}\binom{n}{d}\right)\right)^{2}}{2\left(4+2\log\left(\tfrac{d}{\epsilon}\binom{n}{d}\right)\right)^{2}}\right)$
$\displaystyle=\exp\left(-\frac{\left(2+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)\right)^{2}}{2\left(4+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)\right)}\right).$

Now we take the union bound over every subset of $[n]$ with cardinality $d$
with one specific index chosen. This gives an upper bound for $P_{e}^{n,d}$.

$\displaystyle P_{e}^{n,d}$
$\displaystyle\leq\sum_{\begin{subarray}{c}I\subset[n],|I|=d\\\ i_{1}\in
I\end{subarray}}\operatorname{\textsf{P}}\left(\sum_{k\in[t]}X_{k}<2\right)$
$\displaystyle\leq
d\binom{n}{d}\exp\left(-\frac{\left(2+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)\right)^{2}}{2\left(4+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)\right)}\right)$
$\displaystyle\leq
d\tbinom{n}{d}\exp\left(-1-\log\left(\tfrac{d}{\epsilon}\tbinom{n}{d}\right)\right)$
$\displaystyle\leq
d\tbinom{n}{d}\exp\Big{(}\log\Big{(}\tfrac{\epsilon}{d\tbinom{n}{d}}\Big{)}\Big{)}$
$\displaystyle=\epsilon$

To study the asymptotic behaviour of the construction in Lemma 1, we reduce
our final expression using the following asymptotic approximation of the
binomial coefficient

$\binom{n}{\alpha
n}\approx\frac{1}{\sqrt{2\pi\alpha(1-\alpha)n}}2^{nh(\alpha)},$

where $h(\cdot)$ is the binary entropy function. This gives

$\displaystyle\frac{\binom{n}{\lfloor\frac{n}{d}\rfloor}}{\binom{n-d}{\lfloor\frac{n}{d}\rfloor-1}}\approx\sqrt{\frac{n-d}{n}}2^{dh(1/d)}<\sqrt{\frac{n-d}{n}}ed,$

so that our construction is a $(t,n,d,\epsilon)$-tropical code if

$t>\sqrt{\frac{n-d}{n}}ed\left(4+2\log\left(\frac{d}{\epsilon}\binom{n}{d}\right)\right).$
(5)

Using the bounds $\sqrt{\frac{n-d}{n}}<1$ and
$\binom{n}{d}\leq\left(\frac{en}{d}\right)^{d}$ we can relax the condition in
(5) so that our construction is a $(t,n,d,\epsilon)$-tropical code if

$t>2e\left(d^{2}\log\left(\frac{ne}{d}\right)+d\log\left(\frac{d}{\epsilon}\right)+2d\right).$
(6)

∎

##  V Comparison

Since both probabilistic and deterministic constructions exist for disjunct
matrices, and our result in Theorem 1 provides a construction for a
$(d-1)$-doubly disjunct block design from a $(d-1)$-disjunct block design, we
have three ways to produce $(t,n,d,\epsilon)$-tropical codes for some fixed
$d$:

  1. 1.

First construct a deterministic $(d-1)$-disjunct matrix, then use Theorem 2 to
construct a $(d-1)$-doubly disjunct matrix;

  2. 2.

Construct a $(d-1)$-disjunct matrix probabilistically, then use Theorem 2 to
construct a $(d-1)$-doubly disjunct matrix;

  3. 3.

Construct a $(d-1)$-doubly disjunct random matrix directly using Lemma 1.

Theorem 2 doubles the number of tests used, but does not change the asymptotic
behaviour of the construction used for a $(d-1)$-disjunct matrix. The
deterministic constructions given by [8, 3] have $t=O(d^{2}\log n)$ with
$\epsilon=0$ for specific parameter choices. As such, the first method also
gives $t=O(d^{2}\log n)$ with $\epsilon=0$. For the second method, we can use
our construction detailed in Section IV-B to construct $(d-1)$-disjunct
matrices by replacing the right hand side of Lemma 1 and all subsequent
calculations with $\frac{1+\Delta}{t}$. This gives $(d-1)$-disjunct matrices
with
$t>e\left(2d^{2}\log\left(\frac{ne}{d}\right)+2d\log\left(\frac{d}{\epsilon}\right)+3d\right)$
for error bound $\epsilon$. Then, using Theorem 2, we can get a $(d-1)$-doubly
disjunct matrix with
$t>2e\left(2d^{2}\log\left(\frac{ne}{d}\right)+2d\log\left(\frac{d}{\epsilon}\right)+3d\right)$,
which is greater than the bound in Theorem 3 for the construction of a
$(d-1)$-doubly disjunct matrix directly, as in the third method. From equation
(5) in Theorem 3, the construction given in Section IV-B requires
$t=O\left(d\log\binom{n}{d}+d\log\left(\frac{d}{\epsilon}\right)\right)$. When
$d=o(n)$ is sublinear, then we have $t=O(d^{2}\log n+d\log(d/\epsilon))$. If
$d=\Theta(n)$, then $t=o(d^{2}\log n+d\log(d/\epsilon))$.

In a practical setting, if the situation requires zero error, then we can use
the first method to construct a doubly disjunct matrix with the set of given
parameters. Otherwise, if some error $\epsilon$ is tolerable, then our
construction in Section IV-B may be preferable.

##  VI Conclusion

In this paper we consider a probabilistic construction for a $(d-1)$-doubly
disjunct block design, as well as a deterministic method based on the
existence of a $(d-1)$-disjunct block design. We show that constructing a
$(d-1)$-doubly-disjunct block design probabilistically gives a similar
performance to a fully deterministic construction based on a $(d-1)$-disjunct
block design, with the added benefit of being able to use arbitrary
parameters.

## References

  * [1] H.-P. Wang, R. Gabrys, and A. Vardy, “Tropical group testing,” 2022. [Online]. Available: https://arxiv.org/abs/2201.05440
  * [2] R. Graham and N. Sloane, “Lower bounds for constant weight codes,” _IEEE Transactions on Information Theory_ , vol. 26, no. 1, pp. 37–43, 1980.
  * [3] W. Kautz and R. Singleton, “Nonrandom binary superimposed codes,” _IEEE Transactions on Information Theory_ , vol. 10, no. 4, pp. 363–377, 1964.
  * [4] M. Cheraghchi and V. Nakos, “Combinatorial group testing and sparse recovery schemes with near-optimal decoding time,” _CoRR_ , vol. abs/2006.08420, 2020\. [Online]. Available: https://arxiv.org/abs/2006.08420
  * [5] H. A. Inan, P. Kairouz, M. Wootters, and A. Özgür, “On the optimality of the kautz-singleton construction in probabilistic group testing,” _IEEE Trans. Inf. Theor._ , vol. 65, no. 9, p. 5592–5603, sep 2019. [Online]. Available: https://doi.org/10.1109/TIT.2019.2902397
  * [6] S. Boyd and L. Vandenberghe, _Convex Optimization_. Cambridge University Press, 2004.
  * [7] M. Mitzenmacher and E. Upfal, _Probability and computing: Randomization and probabilistic techniques in algorithms and data analysis_. Cambridge university press, 2017.
  * [8] E. Porat and A. Rothschild, “Explicit non-adaptive combinatorial group testing schemes,” _CoRR_ , vol. abs/0712.3876, 2007. [Online]. Available: http://arxiv.org/abs/0712.3876