Statistical Mechanics of Paraparticles

TL;DR

This paper reviews the foundations of paraparticle statistics, extends previous work, and derives an expression for their heat capacity, aiding potential experimental detection of paraparticles in various systems.

Contribution

It provides a detailed review of paraparticle statistics and derives a new formula for heat capacity, advancing understanding and experimental prospects of paraparticles.

Findings

Derived an expression for paraparticle heat capacity

Extended previous foundational work on paraparticle statistics

Facilitates potential experimental detection of paraparticles

Abstract

Quantum mechanics broadly classifies the particles into two categories: $(1)$ fermions and $(2)$ bosons. Fermions are half-integer spin particles, obeying Pauli's exclusion principle and Fermi-Dirac statistics. Whereas bosons are integer spin particles, not obeying Pauli's exclusion principle, and obeying Bose-Einstein statistics. However, there are two exceptions to this standard case: first, anyons, which exist in 2-dimensional systems, and secondly, paraparticles, which can exist in any dimension. Paraparticles follow their non-trivial parastatistics, obeying their generalised exclusion principle. In this paper, we provide a detailed review of the foundations of paraparticle statistics established in \cite{wang2025particle}. We extend this work further and then derive an important expression for the heat capacity of paraparticles for a specific case, which would provide a handle for the experimental detection of paraparticles in appropriate systems.

Tables

Table 1. Table 1: Examples of R R matrices and their single-mode partition functions z R ​ ( x ) z_{R}(x)

Ex.	1	2	3	4
$R_{cd}^{ab}$	$-{\delta }_{ad}{\delta }_{bc}$	${\delta }_{ad}{\delta }_{bc}{\left(-1\right)}^{{\delta }_{ab}}$	$-{\delta }_{ac}{\delta }_{bd}$	${\lambda }_{ab}{ϵ}_{cd}-{\delta }_{ac}{\delta }_{bd}$
$z_R\left(x\right)$	${\left(1+x\right)}^{ℳ}$	${\left(1+x\right)}^{ℳ}$	$1+ℳx$	$1+ℳx+x^2$

Equations

$$\psi (x_2,x_1,...,x_n)=c\psi (x_1,x_2,...,x_n)$$

$$\psi (x_1,x_2,.,x_n)$$

$$\displaystyle=c^{2}\psi(x_{1},x_{2},.,x_{n})\$$

$${\psi }^I(\{x_i{\}}_{i=1}^n){|}_{x_j\leftrightarrow x_{j+1}}=\sum _J(R{)}_J^I{\psi }^J(\{x_i{\}}_{i=1}^n)$$

$${\hat{\psi }}_{i,a}^{-}{\hat{\psi }}_{j,b}^{+}=\sum _{c,d}{R}_{bd}^{ac}{\hat{\psi }}_{j,c}^{+}{\hat{\psi }}_{i,d}^{-}+{\delta }_{ab}{\delta }_{ij}$$

$${\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,b}^{+}=\sum _{c,d}{R}_{ab}^{cd}{\hat{\psi }}_{j,c}^{+}{\hat{\psi }}_{i,d}^{+}$$

$${\hat{\psi }}_{i,a}^{-}{\hat{\psi }}_{j,b}^{-}=\sum _{c,d}{R}_{dc}^{ba}{\hat{\psi }}_{j,c}^{-}{\hat{\psi }}_{i,d}^{-}$$

$${\hat{e}}_{ij}=\sum _{a=1}^{\mathcal{M}}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}$$

$$[{\hat{e}}_{ij},{\hat{\psi }}_{k,b}^{+}]=\sum _{a=1}^{\mathcal{M}}[{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-},{\hat{\psi }}_{k,b}^{+}]$$

$$=\sum _a[{\hat{\psi }}_{i,a}^{+},{\hat{\psi }}_{k,b}^{+}]{\hat{\psi }}_{j,a}^{-}+{\hat{\psi }}_{i,a}^{+}[{\hat{\psi }}_{j,a}^{-},{\hat{\psi }}_{k,b}^{+}]$$

$$\begin{gathered} =\sum _a({\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{j,a}^{-}-{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-} \\ +{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{+}-{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{j,a}^{-}) \end{gathered}$$

$$=\sum _a\left({\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{+}-{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}\right)$$

$$=\sum _a{\hat{\psi }}_{i,a}^{+}\left(\sum _{c,d}{R}_{bd}^{ac}{\hat{\psi }}_{k,c}^{+}{\hat{\psi }}_{j,d}^{-}+{\delta }_{jk}{\delta }_{ab}\right)-\sum _a{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}$$

$$=\sum _{a,c,d}{R}_{bd}^{ac}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{k,c}^{+}{\hat{\psi }}_{j,d}^{-}+\sum _a{\hat{\psi }}_{i,a}^{+}{\delta }_{jk}{\delta }_{ab}-\sum _a{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}$$

$$=\sum _d\left(\sum _{a,c}{R}_{bd}^{ac}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{k,c}^{+}\right){\hat{\psi }}_{j,d}^{-}+{\delta }_{jk}{\hat{\psi }}_{i,b}^{+}-\sum _a{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}$$

$$=\sum _d{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{i,d}^{+}{\hat{\psi }}_{j,d}^{-}+{\delta }_{jk}{\hat{\psi }}_{i,b}^{+}-\sum _a{\hat{\psi }}_{k,b}^{+}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}$$

$$[{\hat{e}}_{ij},{\hat{\psi }}_{k,b}^{+}]={\delta }_{jk}{\hat{\psi }}_{i,b}^{+}$$

$$[{\hat{e}}_{ij},{\hat{\psi }}_{k,b}^{-}]=\sum _{a=1}^{\mathcal{M}}[{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-},{\hat{\psi }}_{k,b}^{-}]$$

$$=\sum _a[{\hat{\psi }}_{i,a}^{+},{\hat{\psi }}_{k,b}^{-}]{\hat{\psi }}_{j,a}^{-}+{\hat{\psi }}_{i,a}^{+}[{\hat{\psi }}_{j,a}^{-},{\hat{\psi }}_{k,b}^{-}]$$

$$\begin{gathered} =\sum _a({\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{k,b}^{-}{\hat{\psi }}_{j,a}^{-}-{\hat{\psi }}_{k,b}^{-}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-} \\ +{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{-}-{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{k,b}^{-}{\hat{\psi }}_{j,a}^{-}) \end{gathered}$$

$$=\sum _a\left(-{\hat{\psi }}_{k,b}^{-}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}+{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{-}\right)$$

$$=-\sum _a\left(\sum _{c,d}{R}_{ad}^{bc}{\hat{\psi }}_{i,c}^{+}{\hat{\psi }}_{k,d}^{-}+{\delta }_{ba}{\delta }_{ki}\right){\hat{\psi }}_{j,a}^{-}+\sum _a{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{-}$$

$$=-\sum _{a,c,d}{R}_{ad}^{bc}{\hat{\psi }}_{i,c}^{+}{\hat{\psi }}_{k,d}^{-}{\hat{\psi }}_{j,a}^{-}-\sum _a{\hat{\psi }}_{j,a}^{-}{\delta }_{ba}{\delta }_{ki}+\sum _a{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{-}$$

$$=-\sum _c{\hat{\psi }}_{i,c}^{+}\left(\sum _{a,d}{R}_{ad}^{bc}{\hat{\psi }}_{k,d}^{-}{\hat{\psi }}_{j,a}^{-}\right)-{\delta }_{ki}{\hat{\psi }}_{j,b}^{-}+\sum _a{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{-}$$

$$=-\sum _c{\hat{\psi }}_{i,c}^{+}{\hat{\psi }}_{j,c}^{-}{\hat{\psi }}_{k,b}^{-}-{\delta }_{ki}{\hat{\psi }}_{j,b}^{-}+\sum _a{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{j,a}^{-}{\hat{\psi }}_{k,b}^{-}$$

$$[{\hat{e}}_{ij},{\hat{\psi }}_{k,b}^{-}]=-{\delta }_{ki}{\hat{\psi }}_{j,b}^{-}$$

$$[{\hat{e}}_{ij},{\hat{e}}_{kl}]$$

$$=\sum _b[{\hat{e}}_{ij},{\hat{\psi }}_{k,b}^{+}]{\hat{\psi }}_{l,b}^{-}+\sum _b{\hat{\psi }}_{k,b}^{+}[{\hat{e}}_{ij},{\hat{\psi }}_{l,b}^{-}]$$

$$[{\hat{e}}_{ij},{\hat{e}}_{kl}]=\sum _b{\delta }_{jk}{\hat{\psi }}_{i,b}^{+}{\hat{\psi }}_{l,b}^{-}-\sum _b{\hat{\psi }}_{k,b}^{+}{\delta }_{il}{\hat{\psi }}_{j,b}^{-}$$

$$[{\hat{e}}_{ij},{\hat{e}}_{kl}]={\delta }_{jk}{\hat{e}}_{il}-{\delta }_{il}{\hat{e}}_{kj}$$

$${\hat{n}}_i={\hat{e}}_{ii}=\sum _{a=1}^{\mathcal{M}}{\hat{\psi }}_{i,a}^{+}{\hat{\psi }}_{i,a}^{-}$$

Full text

I Introduction

In quantum mechanics, a many-particle state is defined by the wave function $\psi(x_{1},x_{2},...,x_{n})$ which depends on the position coordinates of the particles. When we swap two particles, it does not change the physical state, but it does change the phase of the wavefunction. If we swap a particle at position $x_{1}$ with a particle at position $x_{2}$ we get:

[TABLE]

If we swap the same two particles again, then the system should go back to its original state. Performing a second exchange:

[TABLE]

Since the wavefunction must remain non-zero $\psi\neq 0$, this implies $c^{2}=1$, which gives us two solutions: (1) $c=+1:$ The wavefunction does not change sign when the particles are exchanged. These particles are called bosons. (2) $c=-1:$ The wavefunction acquires a negative sign when particles are exchanged. These particles are called fermions.

Nevertheless, there are two notable exceptions to these quantum exchange statistics: the first is anyons, which can be found solely in two dimensions, and the second is paraparticles that adhere to parastatistics, which can exist in any dimensional space wang2025particle .

For anyons (in 2D) $\psi(x_{2},x_{1},...,x_{n})=e^{i\theta}\psi(x_{1},x_{2},...,x_{n})$, where $e^{i\theta}=$ phase factor and $\theta$ can take the value from [math] to $\pi$.

However, the exchange rules for paraparticles are as follows: the wavefunction for paraparticles is written as $\psi^{I}(\{x_{i}\}_{i=1}^{n})$ where superscript $I$ is a collective index denoting all the flavours (internal quantum numbers) of paraparticles. $x_{1},x_{2},...,x_{n}$ denotes the modes of paraparticles. Modes can be position coordinates, momentum states, etc.

The most fundamental aspect is the exchange rule. When two paraparticles $(x_{j})$ and $(x_{j+1})$ are swapped, the exchange rule is given as:

[TABLE]

$\sum_{J}$: Summation over all the possible flavours (internal quantum states) labelled by $J$ and $|_{x_{j}\leftrightarrow x_{j+1}}$ means swap particles $j$ and $j+1$ (exchange their modes).

$(R)_{J}^{I}$: matrix associated with the exchange operation which reshuffles flavours (internal states) instead of just applying a simple phase change. This is where the paraparticle exchange statistics differ from those of fermions and bosons.

II Commutation Relations for $R$-Matrix

The creation and annihilation operators for paraparticles are defined as:

$\hat{\psi}_{i,a}^{-}$ : Annihilation Operator - Removes a paraparticle in mode $i$ of flavour $a$ (internal quantum number).

$\hat{\psi}_{j,b}^{+}$ : Creation Operator - Adds a paraparticle in mode $j$ of flavour $b$.

For a given $R$-matrix, the paraparticle creation and annihilation operators are defined through the commutation relations wang2025particle :

[TABLE]

[TABLE]

[TABLE]

where $i$ and $j$ are different modes (for example, position coordinates or momentum states) and $a,b,c$ and $d$ are the flavours (for example, spin, etc).

II.1 Contracted Bilinear Operator

A contracted bilinear operator tracks how one type of paraparticle transforms into another, or it allows us to describe how paraparticles move between different modes while keeping track of their flavours. It is defined as:

[TABLE]

This operator obeys Lie algebra rules wang2025particle . The commutative relations of this operator with creation and annihilation operators, as well as with itself, tell us about the transformation of flavours. In the following, we have derived these expressions in detail:

II.1.1 To prove: $[\hat{e}_{ij},\hat{\psi}_{k,b}^{+}]=\delta_{jk}\hat{\psi}_{i,b}^{+}$

Proof:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

The product ${\hat{\psi}_{j,a}^{-}\hat{\psi}_{k,b}^{+}}$ is replaced using equation $(4)$, we get:

[TABLE]

[TABLE]

[TABLE]

$\left(\sum_{a,c}R_{bd}^{ac}\hat{\psi}_{i,a}^{+}\hat{\psi}_{k,c}^{+}\right)$ is replaced using equation (5) and we get:

[TABLE]

Using equation (7), $\sum_{d}\hat{\psi}_{k,b}^{+}\hat{\psi}_{i,d}^{+}\hat{\psi}_{j,d}^{-}$ becomes $\hat{\psi}_{k,b}^{+}\hat{e}_{ij}$ and $\sum_{a}\hat{\psi}_{k,b}^{+}\hat{\psi}_{i,a}^{+}\hat{\psi}_{j,a}^{-}$ becomes $\hat{\psi}_{k,b}^{+}\hat{e}_{ij}$ and the two terms cancel each other and we get:

[TABLE]

II.1.2 To prove: $[\hat{e}_{ij},\hat{\psi}_{k,b}^{-}]=-\delta_{ki}\hat{\psi}_{j,b}^{-}$

Proof:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

The product $\hat{\psi}_{k,b}^{-}\hat{\psi}_{i,a}^{+}$ is replaced using equation $(4)$, we get:

[TABLE]

[TABLE]

[TABLE]

$\left(\sum_{a,d}R_{ad}^{bc}\hat{\psi}_{k,d}^{-}\hat{\psi}_{j,a}^{-}\right)$ is replaced using equation (6) and we get:

[TABLE]

Using equation (7), $\sum_{c}\hat{\psi}_{i,c}^{+}\hat{\psi}_{j,c}^{-}\hat{\psi}_{k,b}^{-}$ becomes $\hat{e}_{ij}\hat{\psi}_{k,b}^{-}$ and $\sum_{a}\hat{\psi}_{i,a}^{+}\hat{\psi}_{j,a}^{-}\hat{\psi}_{k,b}^{-}$ becomes $\hat{e}_{ij}\hat{\psi}_{k,b}^{-}$ and the two terms cancel each other and we get:

[TABLE]

II.1.3 To prove: $[\hat{e}_{ij},\hat{e}_{kl}]=\delta_{jk}\hat{e}_{il}-\delta_{il}\hat{e}_{kj}$

Proof:

[TABLE]

Using equation $(8)$ and $(9)$:

[TABLE]

[TABLE]

II.2 Particle Number Operator

This operator counts the number of particles in a particular mode, summing over flavours. Mathematically:

[TABLE]

The commutation relations for the paraparticle number operator are as follows:

$[\hat{n}_{i},\hat{n}_{j}]=[\hat{e}_{ii},\hat{e}_{jj}]=\delta_{ij}\hat{e}_{ij}-\delta_{ij}\hat{e}_{ji}=0$ 2. 2.

$[\hat{n}_{i},\hat{\psi}_{j,b}^{+}]=[\hat{e}_{ii},\hat{\psi}_{j,b}^{+}]=\delta_{ij}\hat{\psi}_{i,b}^{+}$ 3. 3.

$[\hat{n}_{i},\hat{\psi}_{j,b}^{-}]=[\hat{e}_{ii},\hat{\psi}_{j,b}^{-}]=-\delta_{ij}\hat{\psi}_{i,b}^{-}$ Combining the two: $[\hat{n}_{i},\hat{\psi}_{j,b}^{\pm}]=\pm\delta_{ij}\hat{\psi}_{j,b}^{\pm}$ 4. 4.

Total number operator: $\hat{n}=\sum_{i}\hat{n}_{i}$ , combining this with the second relation we get: $[\hat{n},\hat{\psi}_{j,b}^{+}]=\hat{\psi}_{j,b}^{+}$

II.3 Mathematical Proof of the Exclusion Principle for Paraparticles (Para-exclusion Principle)

For one of the paraparticle statistics cases where there can be only one paraparticle in a given mode we have $\hat{\psi}_{i,a}^{+}\hat{\psi}_{i,b}^{+}=0$ which can be proved by using the specific R statistics (Table 1 in wang2025particle ):

Taking the value of $R_{cd}^{ab}$ from Table 1 (Example 3), $R_{cd}^{ab}=-\delta_{ac}\delta_{bd}$, and putting it in equation (5), we get:

[TABLE]

[TABLE]

Taking $i=j$ (same mode), $\hat{\psi}_{i,a}^{+}\hat{\psi}_{i,b}^{+}=-\hat{\psi}_{i,a}^{+}\hat{\psi}_{i,b}^{+}$ or $2\hat{\psi}_{i,a}^{+}\hat{\psi}_{i,b}^{+}=0$ $\Rightarrow$ $\hat{\psi}_{i,a}^{+}\hat{\psi}_{i,b}^{+}=0$. This means that a single mode (the above formula has used mode $i$) can not be occupied by two paraparticles even if they have different flavours $(a{\neq}b)$.

II.4 State Space

The construction of the state space of the system containing multiple paraparticles is done by applying the creation operator on the vacuum (just like in the case of fermion Fock Space) wang2025particle .

[TABLE]

where, $n$ is the total number of paraparticles, $i$ represents the different modes and $a$ represents the different flavours.

II.5 Unitary Exchange Operator

The initial quantum state refers to a situation where we have two paraparticles located at different positions in space, mathematically: $|0:ia,jb\rangle=\hat{\psi}_{i,a}^{+}\hat{\psi}_{j,b}^{+}|0\rangle$, where, $|0\rangle$ means vacuum state, no particles yet and $\hat{\psi}_{i,a}^{+}$ is defined in starting of section II. The Unitary Exchange Operator $\hat{E}_{ij}$ is defined such that it swaps the position of paraparticles at modes $i$ and $j$ wang2025particle .

[TABLE]

which moves particle $i$ to $j$ and vice versa.

If we apply $\hat{E}_{ij}$ to original state:

[TABLE]

[TABLE]

Using equations $(13)$, and vacuum remains invariant under any operator. Also using equation $(5)$ in next step to change $\hat{\psi}_{j,a}^{+}\hat{\psi}_{i,b}^{+}$ to $\sum_{a^{\prime},b^{\prime}}R^{b^{\prime}a^{\prime}}_{ab}\hat{\psi}_{i,b^{\prime}}^{+}\hat{\psi}_{j,a^{\prime}}^{+}$.

[TABLE]

[TABLE]

Paraparticles swap their positions, and the flavours transform according to the R-matrix. The unitary rotation shows that paraparticles do not simply swap like fermions or bosons, but their states transform non-trivially.

III Exact Solution of a Free Paraparticle

Hamiltonian $\hat{H}$ describes the total energy of the system of free (non-interacting) paraparticles wang2025particle :

[TABLE]

where, $i,j=1,2,...,N$ are the modes (could be lattice sites or momentum states), $a=1,2,...,{\mathscr{M}}$ are the flavours (internal states), $h_{ij}$ is the interaction strength between modes $i$ and $j$ and $h_{ij}=h_{ij}^{*}$ is Hermitian matrix that encodes how paraparticles move between the modes.

We need to diagonalise the Hamiltonian for further computation.

The current form is complicated because paraparticles at different sites are mixed through $h_{ij}$ - the Hamiltonian has off-diagonal terms in $i$, $j$ where modes are coupled. This makes time evolution and observables difficult to compute. We want to diagonalise the Hamiltonian, that is, find a new set of operators where it is decoupled so that separate modes have separate Hamiltonians.

Defining new operators $\hat{\psi}_{i,a}^{-}=\sum_{k=1}^{N}U_{ki}^{*}\tilde{\psi}_{k,a}^{-}$ and $\hat{\psi}_{i,a}^{+}=\sum_{k=1}^{N}U_{ki}\tilde{\psi}_{k,a}^{+}$, where U$$\in{C}^{N\times N}= a unitary matrix that diagonalizes $h$. $UhU^{\dagger}={diag}({\epsilon}_{1},...,{\epsilon}_{N})$ and $UU^{\dagger}=U^{\dagger}U=1$. Matrix h is diagonalised as $h_{kp}=(UhU^{\dagger})_{kp}={\epsilon}_{k}\delta_{kp}$ or $\sum_{ij}U_{ki}h_{ij}U_{pj}^{*}={\epsilon}_{k}\delta_{kp}$, which gives us ${\epsilon}_{k}=\sum_{ij}U_{ki}h_{ij}U_{kj}^{*}$.

After the Transformation, the Hamiltonian becomes

[TABLE]

where, $\epsilon_{k}=$ eigenvalues of $h_{ij}$.

These are the energy levels of the system, each corresponding to a decoupled mode labelled by $k$. Physically, they describe the energy required to occupy the mode $k$. Now $\hat{H}$ has a simple form; each paraparticle has its energy level. The system now behaves like a collection of independent modes. Each mode is characterised by energy $\epsilon$ and can be occupied by a paraparticle of different flavours ${\mathscr{M}}$.

IV Grand Canonical Partition Function

The partition function $Z$ of the whole system of paraparticles is a product of single-mode partition functions $z_{R}(x_{k})$ wang2025particle .

[TABLE]

where, $x_{k}=e^{-\beta(\epsilon_{k}-\mu)}$ is the Boltzmann Factor, $\beta=\frac{1}{k_{b}T}$, $\epsilon_{k}$ is the energy of mode $k$ and $\mu$ is the chemical potential of the system.

The partition function $Z$ encodes all the thermodynamic information. The trace is over the State space of paraparticles. In the case of free paraparticles, the Hamiltonian has been diagonalised into modes, and each mode contributes independently to the partition function $Z$.

Single-mode partition function is given as wang2025particle :

[TABLE]

where $n$ is the number of paraparticles occupying the mode and $d_{n}$** are non-negative integers that define the generalised exclusion statistics for the paraparticles associated with $R$.** For bosons, there is no exclusion, so $d_{n}=1$ for all $n\geq 0$. For fermions, the Pauli exclusion principle allows at most one particle per mode, hence $d_{0}=d_{1}=1$, and $d_{n}=0$ for $n\geq 2$. Mode $k$ defines one particle states. It is different from the wave number $k$, which can accommodate two electrons of opposite spin. From now onwards, $k$ is defined as one particle state; otherwise, it will be explicitly mentioned if it is wavevector $k$, so that confusion should not arise. For paraparticles, the allowed occupation numbers are governed by a specific $R$-matrix, leading to a generalised exclusion rule where $d_{n}$ depends on both $n$ and flavours (internal degrees of freedom). For example for $R$-matrix defined by exapmle 3 in the table (Table 1 below) we have $d_{0}=1,d_{1}={\mathscr{M}}$, and $d_{n}=0$ for $n\geq 2$, where ${\mathscr{M}}$ is the number of flavours of the paraparticles. For $R$-matrix defined by example 4 in the table (Table 1 below), we have $d_{0}=1,d_{1}={\mathscr{M}},d_{2}=1$, and $d_{n}=0$ for $n\geq 3$.

The partition function can be derived as: The Hamiltonian for a system of non-interacting paraparticles is given by: $\hat{H}=\sum_{k}\epsilon_{k}\hat{n}_{k}$, where $\epsilon_{k}$ is the energy of mode $k$, and $\hat{n}_{k}$ is the number operator for mode $k$. For the ease of calculations, we have absorbed the chemical potential $\mu$ into $\epsilon_{k}$. The canonical partition function is $Z=\text{Tr}[e^{-\beta\hat{H}}]=\text{Tr}[e^{-\beta\sum_{k}\epsilon_{k}\hat{n}_{k}}]$. Because the number operators commute $([\hat{n}_{k},\hat{n}_{k^{\prime}}]=0)$, the exponential factors into a product $e^{-\beta\sum_{k}\epsilon_{k}\hat{n}_{k}}=\prod_{k}e^{-\beta\epsilon_{k}\hat{n}_{k}}$. The trace of the product becomes a product of traces over independent modes $Z=\text{Tr}\left[\prod_{k}e^{-\beta\epsilon_{k}\hat{n}_{k}}\right]=\prod_{k}\text{Tr}[e^{-\beta\epsilon_{k}\hat{n}_{k}}]$. The single-mode partition function is given as $z_{k}(x_{k})=\sum_{n=0}^{\infty}d_{n}e^{-\beta\epsilon_{k}n}=\text{Tr}[e^{-\beta\epsilon_{k}\hat{n}_{k}}]$. Therefore, the partition function $Z$ can be written as a product of the single-mode partition function. Hence proved $Z=Tr[e^{-\beta\hat{H}}]=\prod_{k}z_{R}(x_{k})$.

V Average Occupation Number

We derive the average occupation number $\langle\tilde{n}_{k}\rangle$ for a paraparticle defined by $R$ statistics of Example 3 and Example 4 of Table 1 in reference wang2025particle , which we redraw for completeness here.

In Example $1$, the $R$-matrix is given by $-\delta_{ad}\delta_{bc}$. This rule behaves similarly to fermions but for multiple flavours (internal quantum numbers) of fermions. In this system, there are ${\mathscr{M}}$ flavours of fermions, each one is allowed to occupy the mode (one particle state) individually, but no two particles of the same flavour can exist in the same mode. The partition function in this example is given by $(1+x)^{\mathscr{M}}$. In Example $2$, the $R$-matrix is given by $\delta_{ad}\delta_{bc}(-1)^{\delta_{ab}}$. In this case, the behaviour is identical to ${\mathscr{M}}$ flavour fermions if and only if $a=b$. Otherwise, it exhibits non-trivial statistics. However, the same exclusion rule is enforced such that only one particle of any one of the flavour ${\mathscr{M}}$ is permitted in each mode. The partition function in this example is given by $(1+x)^{\mathscr{M}}$. Example $3$ is the most interesting one, the $R$-matrix is given by $-\delta_{ac}\delta_{bd}$. In this case, even when paraparticles are of different flavours, the occupation of the same mode by two paraparticles is prohibited. This represents a stricter condition than fermions. The partition function in this example is given by $(1+{\mathscr{M}}x)$. For each mode, two possibilities are permitted: (1) The mode may remain unoccupied (1 configuration). (2) The mode may be occupied by a single paraparticle of any type (${\mathscr{M}}$ possible configurations). This is a distinct exclusion principle that cannot be reduced to either fermionic or bosonic statistics, representing a nontrivial extension of conventional quantum statistics. In Example $4$ the $R$-matrix is given by $\lambda_{ab}\epsilon_{cd}-\delta_{ac}\delta_{bd}$. The parameters $\lambda$ and $\epsilon$ are defined as special constant matrices. Their primary function is to enable the mathematical formulation of this exotic behaviour. The partition function in this example is given by $(1+{\mathscr{M}}x+x^{2})$. The system exhibits three distinct configurations: (1) A single configuration is permitted for the unoccupied state. (2) ${\mathscr{M}}$ distinct configurations are allowed for single-paraparticle occupation. (3) One specially constrained configuration is enabled for two-paraparticle occupation. Therefore the partition function is $(1+{\mathscr{M}}x+x^{2})$. The system is allowed to have two paraparticles in the same mode, but only in one very specific way, determined by the R-matrix. This behaviour is fundamentally different from both: (1) Bosonic systems, where particles can freely occupy the same mode in multiple indistinguishable ways. (2) Fermionic systems, where no two particles can occupy the same mode at all due to the Pauli exclusion principle. For systems with ${\mathscr{M}}\geq 3$, the $R$-matrix becomes non-unitary, and the corresponding dynamics may no longer be Hermitian. This indicates a deviation from conventional quantum statisticswang2025particle .

The average occupation number is given by:

[TABLE]

$\frac{d}{dx_{k}}z_{R}(x_{k})=\frac{d}{dx_{k}}\left(\sum_{n=0}^{n_{\max}}d_{n}x_{k}^{n}\right)=\sum_{n=0}^{n_{\max}}nd_{n}x_{k}^{n-1}$. Multiplying both sides by $x_{k}$: $x_{k}\frac{d}{dx_{k}}z_{R}(x_{k})=\sum_{n=0}^{n_{\max}}nd_{n}x_{k}^{n}$ which gives us:

[TABLE]

V.0.1 Case 1: Fermions

Partition function of fermions is $(1+x)$, so the average occupation number is given by:

[TABLE]

V.0.2 Case 2: Bosons

Partition function of fermions is ${(1-x)}^{-1}$, so the average occupation number is given by:

[TABLE]

V.0.3 Case 3: Taking Example 2 from Table 1

Partition function of fermions is ${(1+x)}^{\mathscr{M}}$, so the average occupation number is given by:

[TABLE]

[TABLE]

V.0.4 Case 4: Taking Example 3 from Table 1 (maximum 1 paraparticle per mode)

One of the most interesting cases is that of Example $3$. The mode configurations are as follows: $d_{0}=1$ (vacuum state, no paraparticle), $d_{1}={\mathscr{M}}$ (one paraparticle which can have any one of the ${\mathscr{M}}$ flavours), and $d_{n}=0$ for $n\geq 2$ (This means there can’t be two paraparticles in one mode even if their flavours are different). For this specific case the partition function is given by $z_{R}(x)=1+{\mathscr{M}}x$, where $x=e^{-\beta(\epsilon_{k}-\mu)}$. The average occupation number is given by:

[TABLE]

[TABLE]

This differs from the case of fermions by the factor $\frac{1}{{\mathscr{M}}}$ multiplied to $e^{\beta(\epsilon_{k}-\mu)}$.

V.0.5 Case 5: Taking Example 4 from Table 1 (maximum 2 paraparticles per mode)

Partition function is given by: $z_{R}(x)=1+{\mathscr{M}}x+x^{2}$

[TABLE]

VI Thermodynamics of the Paraparticles of Example 3

The partition function of mode $k$ is:

[TABLE]

The Thermodynamic Potential (analogue of the Gibbs free energy) for $k^{th}$ mode is:

[TABLE]

Therefore,

[TABLE]

[TABLE]

Thermodynamic Potential of the whole gas:

[TABLE]

Entropy of paraparticles:

The entropy is given by pathria :

[TABLE]

[TABLE]

Let us redefine $\langle\tilde{n}_{k}\rangle=n_{k}$ for simplicity. Now,

[TABLE]

or

[TABLE]

[TABLE]

we have

[TABLE]

Then

[TABLE]

[TABLE]

Heat capacity at constant volume and chemical potential $(\mu)$ is given by:

[TABLE]

[TABLE]

From Equation (31):

[TABLE]

[TABLE]

We had

[TABLE]

[TABLE]

[TABLE]

From equation (33):

[TABLE]

[TABLE]

Therefore

[TABLE]

[TABLE]

Introducing the paraparticle density of states (refer to Appendix A):

[TABLE]

Let us study the behaviour in the low temperature limit.

$k_{B}T\ll\mu$ the factor $n(\epsilon){(1-n(\epsilon))}$ acts like a delta function.

[TABLE]

Let us evaluate the integral and substitute:

[TABLE]

[TABLE]

[TABLE]

Consider the integral $I_{1}$

[TABLE]

[TABLE]

[TABLE]

[TABLE]

Let

[TABLE]

[TABLE]

The expression in each parenthesis is the sum of a geometric progression. In total, we have,

[TABLE]

Going back to $x$:

[TABLE]

Therefore,

[TABLE]

[TABLE]

We use the standard integral

[TABLE]

Therefore,

[TABLE]

[TABLE]

Here $Li_{2}$ is the dilogarithmic function. In general $Li_{s}(z)=\sum_{k=1}^{\infty}\frac{z^{k}}{k^{s}},\quad|z|<1$.

Refer back to equation $(40)$ and define

[TABLE]

[TABLE]

[TABLE]

Let

[TABLE]

[TABLE]

[TABLE]

Going back to x:

[TABLE]

It is exactly the same expression as in equation $(42)$ except ${\mathscr{M}}$ going to $\frac{1}{{\mathscr{M}}}$. Therefore, the final expression of the heat capacity is as follows:

[TABLE]

Let us check whether this expression goes to the standard expression for electrons for ${\mathscr{M}}=1$, we have:

[TABLE]

Therefore,

[TABLE]

We notice that in the case of ${\mathscr{M}}=1$ we do get back the standard case of electrons, as it should. Let us compare the ratio of the heat capacity of paraparticles to that of the electrons:

[TABLE]

Here $\zeta$ is the Riemann Zeta function. The above equation is our final result. It provides a handle on the thermodynamic detection of the presence of paraparticles, by just doing the careful heat capacity measurements in appropriate systemswang2025particle . Similar expression can be derived for Example $4$ in Table $1$, but it turns out to be very cumbersome.

VII Conclusion

We have reviewed the basic foundations regarding paraparticles laid in Ref.wang2025particle . We extended that work by computing various thermodynamical quantities for paraparticles including entropy. We also computed heat capacity of paraparticles for the most interesting example, that is, Example $3$ in Table $1$. We compared that expression with that of the standard expression for fermions. We noticed that there is an excess heat capacity in the case of paraparticles which is proportional to the square of the logarithm of the flavour index ${\mathscr{M}}$ of the paraparticles. We propose that our expression can provide a foundation for the experimental detection of paraparticles.

VIII Acknowledgement

The authors are thankful to Zhiyuan Wang for correspondence and very useful comments.

Appendix A Paraparticle Density of States

Let us define paraparticle density of states in 3D as $\rho(\epsilon)$. Energy of paraparticle in mode $k$ is given by the dispersion relation $\epsilon_{k}=\frac{\hbar_{k}}{2m}$wang2025personal . The number of states in a spherical shell of radius $k$ and thickness $dk$ is given by:

[TABLE]

where V represents the total volume of the system and ${\mathscr{G}}$ is the degeneracy that defines how many paraparticles are allowed to occupy a mode (for example 3 ${\mathscr{G}}=1$ ). The density of states in $k$-space is given by $\rho(k)=\frac{dN}{dk}$. On substituting, it gives:

[TABLE]

The total number of particles at $T=0$ can then be computed as:

[TABLE]

Here $k_{GW}$ is a general wavevector which generalises the definition of the Fermi wavevector. Solving for the general wavevector $k_{GW}$ gives:

[TABLE]

To compute the density of states in terms of energy, use the relation wang2025personal :

[TABLE]

Thus, $\rho(k)\,dk={\mathscr{G}}\frac{V}{2\pi^{2}}\cdot k^{2}\,dk=\frac{{\mathscr{G}}V}{2\pi^{2}}\left(\frac{2m\epsilon}{\hbar^{2}}\right)\left(\frac{m}{\hbar^{2}}\right)\left(\frac{\hbar^{2}}{2m\epsilon}\right)^{1/2}d\epsilon$. Simplifying: $\rho(\epsilon)={\mathscr{G}}\frac{Vm}{2\pi^{2}\hbar^{2}}\left(\frac{2m}{\hbar^{2}}\right)^{1/2}\epsilon^{1/2}$

Which can also be written as

[TABLE]

here,

[TABLE]

where $m$ is the mass of the paraparticle.