# Precise simulation of single-hole spin control in quantum dots

YuanDong Wang  Department of Physics, Renmin University of China, Beijing
100872, China  JingHan Ni  Department of Physics, Renmin University of China,
Beijing 100872, China  JianHua Wei<EMAIL_ADDRESS>Department of Physics,
Renmin University of China, Beijing 100872, China

###### Abstract

The precise simulation of the initialization, control, and read-out of a
single-hole spin qubit is investigated by accurately solving the extended
Anderson impurity model in the real time domain with the hierarchical
equations of motion approach. The initialization is realized by the ionization
of an exciton with high fidelity. Then, a SU(2) control is achieved via the
combination of Larmor procession of the hole spin in Voigt geometry magnetic
field and rotation about the optical axis with a geometric phase induced by a
picosecond laser pulse. Finally, the read-out of the qubit is implemented
through photocurrent recording. Our theory not only reproduces the recent
experimental results with one set of internal parameters, but also predicts a
maximal fidelity by adjusting the dot-electrode coupling strength.

###### pacs:

72.15.Qm,73.63.Kv,73.63.-b

The trapped single spin in quantum dots (QDs) is a promising qubit which can
be optically controlled within picosecond scale Brunner  _et al._ (2009);
Press  _et al._ (2008); De Greve  _et al._ (2011). Therefore, it is a good
candidate for integrated circuit in quantum information processing(QIP) with
the mature processing technology of semiconductors. In literatures, many
achievements have been made in electron spin qubit in QDs Press  _et al._
(2008); Vamivakas  _et al._ (2010); Delteil  _et al._ (2014). However, the
non-Markovian hyperfine interaction induces decoherence and drops the fidelity
of electron spin control Vink  _et al._ (2009). The valence band holes possess
p-type wave function that leaves small residual dipolar interactionFischer
_et al._ (2008), which highly suppresses contact hyperfine interaction and
induces longer lifetimes of hole spin than that of electron spin Fallahi,
Yılmaz, and Imamoğlu (2010); Chekhovich  _et al._ (2011). Recently, many
experiments on single-hole spin qubit have been preformed, including
initialization, coherent control and read-out Ardelt  _et al._ (2015); De
Greve  _et al._ (2011); Gerardot  _et al._ (2008); Godden  _et al._ (2012).

For the realistic application of hole spin qubit, the high fidelity during the
initializing process is the key requirement in any QIP protocol. Among the
present methods of initialization, the ionization of an exciton has distinct
advantages, which achieved the fidelity of $98.5\%$ Ardelt  _et al._ (2015).
By reducing fine-structure splitting or applying a magnetic field parallel to
growth direction, the fidelity high to $99\%$ is reachable Godden  _et al._
(2010); Brash  _et al._ (2015). Moreover, it is fast (in $p$s) enough to meet
the requirement that the initialization time should $10^{-4}$ order smaller
than decoherence time. As a comparison, the optical pumping Press  _et al._
(2008), one of other methods, can only reach the fidelity of $95\%$ in the
time scale of $\mu$s due to the several loops to prepare a polarized spin
state Heiss  _et al._ (2008).

Whereas above experimental investigations have been actively performed,
theoretical works are not sufficient so far, and the rate equation is the
commonly used method to simulate the hole spin manipulation Breuer and
Petruccione (2002). We comment that the rate equation is not accurate enough
basing on the following two facts. Firstly, in the QDs-based hole spin qubit,
the QD directly couple to metal leads (electrodes) which inevitably impacts on
the qubit. Therefore, what we deal with is a typical quantum open system with
infinity degree of freedoms of the total density matrix, while the rate
equation only concerns the diagonal terms of the reduced density matrix and
treats the dot-electrode couplings by low-order perturbation schemes.
Secondly, the hole spin qubit system is a typical strongly correlated open
system with degree of freedoms of the electron-electron ($e-e$), hole-hole
($h-h$) and electron-hole ($e-h$) interactions, while the rate equation either
neglects this important interaction or treats it in the single electron level.

Obviously, for the theoretical study on hole spin qubit, a non-perturbative
approach is highly required to deal with the basic quantum model involving
different Coulumb interactions. The hierarchical equations of motion (HEOM)
approach we used in the present work can meet this requirement, which
nonperturbatively resolves the combined effects of dot-electrode dissipation,
Coulumb interactions, and non-Markovian memory Jin, Zheng, and Yan (2008); Li
_et al._ (2012); Ye  _et al._ (2016). In this paper, we start from the
extended Anderson impurity model to describe the hole spin qubit, fully
considering the Coulumb interactions interaction and the dot-electrode
couplings. We deal with this quantum model non-perturbatively in the real time
domain to precisely simulate the single-hole spin manipulation.

In what follows, via HEOM approach, the whole process of QIP including
initialization, coherent control and read-out will be precisely simulated. As
will be demonstrated, our theory not only reproduces the recent experimental
results in Ref. Godden  _et al._ (2012) with one set of internal parameters,
but also predicts a maximal fidelity by adjusting the dot-electrode coupling
strength. The complete process of single hole-spin initialization, coherent
control and read-out process are sketched in Fig. 1(a) and (b). In order to
initialize a single-hole spin , a $\sigma_{+}$ resonant circularly polarized
pulse with pulse area of $\pi$ is used which creates a $e-h$ pair that driving
the ground state $\left|cgs\right>$ into a neutral exciton state
$\left|X^{0}_{\uparrow\Downarrow}\right>$. Due to the much larger effective
mass of holes, the hybridization strength of them is much smaller than that of
electrons, i.e. $\Delta_{H}\ll\Delta_{E}$. As a consequence, electrons in
conduction level tunnels into electrode in a $2\sim 3$ orders of time
magnitude faster than holes, which turns
$\left|X^{0}_{\uparrow\Downarrow}\right>$ into single-hole spin state
$\left|\Downarrow\right>$ quickly. Note that the operation is confined to
photocurrent region, we omit the radiative recombination of $e-h$ pair. Then,
the applied in-plane magnetic field $B_{x}$ drives the single-hole spin to
precess along $x$ axis, which preforms a U(1) operation. In order to realize a
SU(2) operation in Bloch sphere, a geometric phase approach is adopted as
proposed in Ref.Economou and Reinecke (2007). That with a geometric pulse of a
sech envelope, the hole spin undergoes a cycle from single hole state to
positive trion state and back to single hole state acquired an rotation angle
$\phi_{z}$. After the control sequence, a detection circular polarized pulse
$\sigma_{+}$ is applied, which partially excites the single-hole spin state to
positive trion state, accompanied with a photocurrent proportional to the
spin-down component of the hole spin state, for which the qubit read-out is
achieved via photocurrent detection technique Mar  _et al._ (2013).

Figure 1: (Color online). (a) Schematic illustration of the single-hole spin
control process. (b) Initialization, coherent control and read-out of a
single-hole spin, where the full arrows represent the optical excitation
process and the dash arrows indicate the transitions due to tunneling. Figure
2: (Color online).Initialization process, with a $\sigma_{+}$ circular
polarized pulse tuned on resonance pules applied, for which the pules area is
$\pi$ and $0\leq t\leq 10$ ps duration. (a) The time evolution of occupation
numbers of $c$-level and $v$-level. (b) The real time photocurrent, with the
maximum appeared at the end of the pulse ($t\sim 10$ ps). (c) Fidelity as a
function of time $t$ which contains the excitation-tunneling process, with
increasing the hybridization strength of hole. (d) Fidelity $F$ as a function
of hybridization strength of $\Delta_{E}$, at which the maximum fidelity
occurs at $\Delta_{E}\simeq 0.5$. The inset shows the half of Rabi cycle $T/2$
as a function of $\Delta_{E}$. (e) Density of spectral function of electrons
with different electron-electrode hybridization strength.

By reference of experimental structures, our single-hole spin qubit consists
of a QD connecting to two electrodes, which can be described by an extended
Anderson impurity Hamiltonian with $e-e$ , $h-h$ and $e-h$ interactions
considered. The total Hamiltonian is written as

$H=H_{c}+H_{v}+H_{c-v}+H_{opt}+H_{res}+H_{dot-res},$ (1)

where $H_{c}$ and $H_{v}$ describe the conduction and valence level with $e-e$
and $h-h$ interactions respectively

$\displaystyle
H_{c}=\sum_{\mu}\epsilon_{c}\hat{n}_{c\mu}+U_{c}\hat{n}_{c\uparrow}\hat{n}_{c\downarrow},$
(2) $\displaystyle
H_{v}=\sum_{\mu}\epsilon_{v}\hat{n}_{v\mu}+U_{v}(1-\hat{n}_{v\uparrow})(1-\hat{n}_{v\downarrow}).$
(3)

In above equations, $\hat{n}_{c\mu}=\hat{a}_{c\mu}^{\dagger}\hat{a}_{c\mu}$
and $\hat{n}_{c\mu}=\hat{a}_{c\mu}^{\dagger}\hat{a}_{c\mu}$, where
$\hat{a}_{c\mu}(\hat{a}_{c\mu}^{\dagger})$ annihilates(creates) a conduction
level of spin $\mu$, and similar to the valence level. $U_{c}(U_{v})$ is the
Coulomb repulsion energies if the $c$-level ($v$-level) is double (zero)
occupied. The term

$H_{c-v}=-\sum_{\mu,\mu^{{}^{\prime}}}U_{exc}\hat{n}_{c\mu}(1-\hat{n}_{v\mu^{{}^{\prime}}})$
(4)

accounts for the Coulomb attraction energies between the $e-h$ pair. $H_{opt}$
denotes the interaction of the control field on QD, whose explicit expression
will be specified later. The electrodes are modeled by non-interaction
electrons

$H_{res}=\sum_{\alpha k\mu}(\epsilon_{\alpha k}+\mu_{\alpha})\hat{d}_{\alpha
k\mu}^{\dagger}\hat{d}_{\alpha k\mu},$ (5)

where $\hat{d}_{\alpha k\mu}(\hat{d}_{\alpha k\mu}^{\dagger})$ denotes the
creation(annihilation) operator of electron in the specified
$\alpha$-electrode spin-orbital state $\left|k\right>$ of energy
$\epsilon_{\alpha k}$. Nonequilibrium chemical potential $\mu_{\alpha}$ with
$\alpha=L,R$ will arise in the presence of bias of voltage. The zero-energy
point is set to be at the equilibrium chemical potential
$\mu_{\alpha}^{eq}=0$. The coupling between the dot and the electrode is
described by

$H_{dot-res}=\sum_{\alpha k\mu}(t_{c\alpha
k}\hat{a}_{c\mu}^{\dagger}\hat{d}_{\alpha k\mu}+t_{v\alpha
k}\hat{a}_{v\mu}^{\dagger}\hat{d}_{\alpha k\mu}+H.c.),$ (6)

It should be noted that due to the large effective mass of holes, transfer
matrix element for conduction level $t_{v\alpha k}$ is assumed to be much
smaller than that of valence level $t_{c\alpha k}$. In the HEOM theory, the
influence of electron reservoirs on the dot acts through the hybridization
functions with a Lorentzian form
$\Delta_{\mu\nu}(\omega)\equiv\sum_{\alpha}\Delta_{\alpha\mu\nu}(\omega)=\pi\sum_{\alpha
k}t_{\alpha\mu k}t_{\alpha\nu k}^{*}\delta(\omega-\epsilon_{\alpha k})=\Delta
W^{2}/[(\omega-\mu_{\alpha})^{2}+W^{2}]$, where $W$ is the bandwidth and
$\mu_{\alpha}$ is the chemical potential of lead $\alpha$. The details of HEOM
formalism has been developed in Refs. Jin, Zheng, and Yan (2008); Li  _et al._
(2012), and the final HEOM can be cast into a compact form of

$\displaystyle\dot{\rho}^{(n)}_{j_{1}\cdots j_{n}}=$
$\displaystyle-\Big{(}i{\cal
L}+\sum_{r=1}^{n}\gamma_{j_{r}}\Big{)}\rho^{(n)}_{j_{1}\cdots
j_{n}}-i\sum_{j}\\!{\cal A}_{\bar{j}}\,\rho^{(n+1)}_{j_{1}\cdots j_{n}j}$
$\displaystyle-i\sum_{r=1}^{n}(-)^{n-r}\,{\cal
C}_{j_{r}}\,\rho^{(n-1)}_{j_{1}\cdots j_{r-1}j_{r+1}\cdots j_{n}},$ (7)

where $\rho_{0}(t)=\rho(t)=tr_{res}\rho_{total}(t)$ is the reduced density
matrix and ${\rho_{j_{1}...j_{n}}(t)}$ are auxiliary density matrices at the
$n$th-tier. Any observable $\hat{O}$ of the dot system can be calculated in
form of $\bar{O}=tr(\rho_{0}\hat{O})$. The transient current through the
electrode $\alpha$ is determined exclusively by the first-tier auxiliary
density operators as,

$I_{\alpha}(t)=e\frac{i}{\hbar^{2}}\sum_{i\mu}tr_{s}\\{\rho_{\alpha\mu}^{\dagger}(t)\hat{a}_{i\mu}-\hat{a}_{i\mu}^{\dagger}\rho_{\alpha\mu}^{-}(t)\\},$
(8)

where the index $i$ sums from $c$ to $v$ that counts the contributions both of
$c$\- and $v$-level. To simulate experiments, we choose the parameters in
$H_{sys}$ having the same energy level structure with experiments, as
schematically shown in Fig. 1(a), with $\epsilon_{c}=2$ (the energy unit is
all set to meV in rest of this paper), $\epsilon_{v}=-2$, $U_{c}=U_{v}=2$,
$U_{exc}=1$, and a reverse bias $V=\mu_{L}-\mu_{R}=0.2$ is applied.

For initialization, a resonant pules with $\pi$ area of 10 ps is applied, with
the Hamiltonian $H_{opt}=\Omega(e^{i\omega
t}c_{c\uparrow}^{\dagger}c_{v\downarrow}+h.c.)$, where $\omega=3$,
$\Omega=0.1$, $\Delta_{E}=0.05$ and $\Delta_{H}=0.0003$. The time evolution of
occupation numbers of $c$\- and $v$-level is presented in Fig. 2(a). As shown
in the figure, the maximum of the electron numbers in $c$-level is around 0.62
rather than unity, which is induced by the fast tunneling electron tunneling
of the neutral exciton $\left|X_{0}\right>$. It will lead to the intensity
damping of the Rabi oscillation with increasing electron-electrodes
hybridization strength $\Delta_{E}$, as a signal of tunneling-induced
dephasing Ardelt  _et al._ (2015).

Figure 3:  (a) The real time photocurrent $I_{\pm}$ in the presence of
magnetic field $B_{x}=0.2$. Initialization pulse ends at $t=10$ ps . With a
co-circular excitation detection pulse applied at the $\tau_{d}=\pi/f_{L}$ ,
no detection photocurrent produced, while with a cross-circular pulse,
detection current is produced. (b) Comparison of $z$ component of the single-
hole spin calculated by $\bar{S}_{z}$ and $\tilde{S}_{z}$ . (c) The schematic
diagram of SU(2) control of a hole spin. the hole spin is initialization to
spin up state and preforms Lamor procession about $x$ axis. The hole spin
points along to $+y$ direction on arrival of the geometric-phase control
pulse. (d) Photocurrent difference $\Delta_{I}=I_{+}-I_{-}$ of $\sigma_{\pm}$
detection pulse as function of precession time. $z$ axis rotation control is
achieved via geometric-phase pulse, where the rotation angle is reflected via
the photocurrent oscillation amplitude.

After the pulse applied, it can be seen that the occupation decay is
exponential for electrons while approximately linear for holes at $t>10$ ps,
due to their different hybridization strength to electrodes. The
initialization is achieved at $t\sim 70$ ps, where the electrons in $c$-level
has totally escaped into electrodes ($N_{c}\sim 0$) and a single hole in
$v$-level has been left ($N_{v}\sim 1$). Fig. 2(b) depicts the photocurrent
along initialization process. The peak of the photocurrent is shown at the end
of the pulse ($t\sim 10$ ps), where the electrons in $c$-level has maximally
accumulated. The charge-photocurrent is precisely conserved, which is
responsible for the photocurrent several orders of magnitude larger than the
experimental data.

In experiments, in order to achieve high fidelity initialization, an AlGaAs
barrier is applied to tailor the tunneling rate of electrons and holes, which
effectively changes the conduction and valence level-electrodes hybridization
strength. To check this effect theoretically, we investigate the dynamics of
fidelity during the excitation-tunneling process at different $\Delta_{E}$ and
$\Delta_{H}$. The fidelity is defined via
$F=\overline{\left<\psi_{in}\right|\hat{U}\hat{\rho}_{out}\hat{U}\left|\psi_{in}\right>}$
Poyatos, Cirac, and Zoller (1997), which measures the distance between the
real evolution $U$ and the target evolution $U_{t}$ to a given initial state
$\left|\psi_{in}\right>$. Here the initial state is set to $\left|cgs\right>$,
the target state is the single hole state $\left|\Downarrow\right>$, which is
arrived with two-step (see Fig. 1(b)) process: a) optical excite the crystal
ground state $\left|cgs\right>$ to neutral exciton $\left|X_{0}\right>$, b) by
fast electron tunneling of the neutral exciton $\left|X_{0}\right>$ to single
hole spin $\left|\Downarrow\right>$. For ideal spin qubit storage with high
fidelity, ultrafast $c$-level electron lifetime and long $v$-level hole
storage against filling from electrodes are required. Fig. 2(c) shows the time
evolution of fidelity with different hole coupling $\Delta_{H}$. Noting that
the electrons in electrodes fill the hole through the whole process, even at
the beginning of the optical excitation. With decreasing hole coupling
$\Delta_{H}$, storage time for hole increases accordingly. At
$\Delta_{H}=0.0003$, as much as $93\%$ of fidelity is observed. This fidelity
dynamics can be observed via the photocurrent amplitude of the trion
transition $X_{+}$ with a probe pulse, which reflects the population of single
hole spin $\left|\Downarrow\right>$.

Ideally, one expect a high initialization fidelity with a ultrafast electron
tunneling. However, when the electron hybridization strength is comparable to
the frequency of the laser which is tuned on resonance with the neutral
excitation, it will inevitably bring damage to fidelity. One of these is the
tunneling-induced dephasing. Ardelt $et~{}al.$ extract this phenomenon in hole
spin initialization with low temperature Ardelt  _et al._ (2015). We comment
that compare to rate equation method, the tunneling-induced dephasing is only
the consequence of the increasing hybridization strength of electrons
$\Delta_{E}$ in HEOM calculations, but not as a parameter.

The max fidelity appearing in real time scanned with $\Delta_{E}$ is shown in
Fig. 2(d). With increasing $\Delta_{E}$, rather than monotonically increase,
the fidelity experiences a maximum with $92.1\%$ at $\Delta_{E}=0.052$. It
should be noted that the couping $\Delta_{E}$ can induce a energy level shift
(the details analysis will be done in Fig. 2(e)) , correspondingly, in order
to maximize the fidelity, the frequency $\omega$ of Rabi oscillation should be
adjusted as well to satisfy resonance condition. The half of Rabi cycle as a
function of $\Delta_{E}$ is shown in the inset of Fig. 2(d). From the
beginning with $\Delta_{E}=0.02$, as $\Delta_{E}$ increases, fidelity starts
to grows fast at first. This growth is dominated by the shorter electron
tunneling time during which the filling of hole is relatively suppressed.
However, continue increase $\Delta_{E}$ that exceeding
$\Delta_{E}/\Delta_{H}\simeq 173$ will bring damage to fidelity, the fidelity
drops to lower than $90\%$ with $\Delta_{E}$ exceeding $0.12$. In order to
explore above impact of electron-electrodes hybridization strength
$\Delta_{E}$ on fidelity, we calculate the single particle spectral function
$A(\omega)$ of $c$-level electron with $\Delta_{E}=0.02$, $0.05$ and $0.1$,
presented in Fig. 2(e). Clearly, due to larger coupling $\Delta_{E}$, the
spectral function $A(\omega)$ experiences linewidth broadening. The broadening
effect can make the single particle level invisible to optical excitation. The
ground state $\left|cgs\right>$ only can be partly pumped to neutral exciton
state $\left|X_{0}\right>$ with a adjacent area around zero detuning. Noting
that the level broadening plays a role as tunneling-induced damping, as the
experiment mentioned above, which can be observed via the broadening of
photucurrent absorption spectral Ardelt  _et al._ (2015). The second feature
of spectral function is the central position shift. At a small couping
$\Delta_{E}=0.02$, the central position is about $\omega\simeq 2.0$, that is
exactly the single particle excitation energy with a $c$-level electron
tunneling event. When $\Delta_{E}$ increases to $0.1$, the central position
moves to $2.08$, meanwhile the fidelity drops to $90.6\%$.

A U(1) rotation of the initialized hole spin state is achieved by applying a
magnetic field $B_{x}$ perpendicular to growth direction to form the Voigt
geometry. The down-spin state $\left|\Downarrow\right>$ is a superposition of
the eigenstates, which will preform Larmor precession about $B_{x}$ with
frequency $f_{L}$ determined by the hole Zeeman splitting. To detect the
single-hole spin, a co-(cross)-circularly polarized detection pulse with $\pi$
pulse area is used to excites the single hole to positive trion state
$\left|X_{+}\right>$ after a time delay $\tau_{d}$. The resonance frequency of
the detection pulse is positively detuned with $U_{exc}$ compared to the
initialization pulse, due to the additional Coulomb interaction between the
$e-h$ pair. To demonstrate the photocurrent detection of a hole spin in our
simulation, we depict the time evolution of photocurrent, the detection pulse
is applied at time delay $\tau_{d}=\pi/f_{L}$, for which the initialized spin
down hole $\left|\Downarrow\right>$ precessed to spin up state
$\left|\Uparrow\right>$. At this time, the excitation of co-circular detection
pulse is completely suppressed nevertheless the cross-circular pulse is
optical active, which is shown in Fig. 3(a). where $\sigma_{+}$ and
$\sigma_{-}$ denoting the co-circular and cross-circular excitation
respectively.

In experiments, the polarization of the hole spin can be read out via
photocurrent using $\tilde{S}_{z}=\frac{I_{-}-I_{+}}{I_{-}+I_{+}}$Godden  _et
al._ (2012), where $I_{\pm}$ is the amplitude of photocurrent peak for the
detection pulse $\sigma_{\pm}$. Theoretically, the exact $z$ component of hole
spin can be directly obtained in form of $\bar{S}_{z}$=tr$(\rho\hat{S}_{z})$,
which can be used to examine the accuracy of the photocurrent read-out
process. Fig. 3(b) shows the comparison between $\tilde{S}_{z}$ and
$\bar{S}_{z}$, with hole hybridization strength $\Delta_{H}=0.0003$, where the
overall consistence between $\tilde{S}_{z}$ and $\bar{S}_{z}$ can be seen. Our
results of $\tilde{S}_{z}$ is in agreement with the experiment data in
Ref.Godden  _et al._ (2012), therefore, the accuracy of single-hole spin
photocurrent read-out technique is verified by Fig. 3(b). When the Rabi
frequency $f_{R}$ and Zeeman splitting energy of the $B_{x}$ magnetic field is
comparable, there exist a small phase difference as shown in the figure, which
indicates that the time evolution of $\tilde{S}_{z}$ is slightly left behind
$\bar{S}_{z}$. The phase difference results from the delay of recording
photocurrent according the definition of $\tilde{S}_{z}$. As shown in Fig.
3(a), the photocurrent reaches its peak value at the end of the pulse, when
the hole spin has undergone a precession already. For the parameters used
here, Rabi frequency $f_{R}=0.1$ and Zeeman energy $B_{x}=0.2$, the phase
delay is of same magnitude as half of Rabi cycle $T_{R}/2\sim 5$ps, as
presented in Fig. 3(b).

In order to implement a SU(2) control of the hole spin, a control pulse is
needed to rotate the spin along second rotation axis. Here we use geometric
phase approach, as proposed in theories Economou  _et al._ (2006); Economou
and Reinecke (2007), and then successfully realized in experiments Godden  _et
al._ (2012); Kim  _et al._ (2010); Greilich  _et al._ (2009). The control
pulse is shaped with a hyperbolic secant envelope, $H_{opt}=\Omega{\rm
sech}(\sigma t)(e^{i\omega t}c_{c\uparrow}^{\dagger}c_{v\downarrow}+h.c.)$
with fixed $\sigma=\Omega=0.2$ to guarantee no population transferring to
trion state after application of the pulse, where $\Omega$ is the Rabi
frequency and $\sigma$ denotes the bandwidth of the pulse. The single-hole
spin state acquires a phase factor
$\phi_{z}=\arctan(\frac{2\sigma\Delta}{\Delta^{2}-\sigma^{2}})$ about $z$-axis
via varying detuning $\Delta$ from the resonance between the single hole spin
state and the positive trion state Economou  _et al._ (2006), as schematically
illustrated in Fig. 1(b).

To be concrete, the geometric-phase control pulse is applied after a time
delay $\tau=27$ ps that the hole spin pointing along +$y$ axis that the pulse
has a maximum effect on the hole spin since the rotation radius is equal to
the radius of Bloch sphere, as shown in Fig. 3(c). Carries with the rotation
angle $\phi_{z}$, the hole spin continues to precess along $x$ axis under the
$B_{x}$ magnetic field, for which the magnitude of geometric phase $\phi_{z}$
determines the rotation radius along $x$ axis. Since the magnitude of the
detection photocurrent is proportional to the hole spin projection on $z$
axis, for which one can use the photocurrent difference $\Delta I=I_{-}-I_{+}$
between $\sigma_{\pm}$ detection pulse to pick up the information of the
geometric rotation angle $\phi_{z}$. The hole spin procession with different
detuning $\Delta$ is shown in Fig. 3(d), that is consistent with the
experiments in Ref.Godden  _et al._ (2012), where the $\sigma\pm$ detection
pulses are scanned through the procession. The detunning are set as
$\Delta=0,\sigma,16$ respectively. For the detunning $\Delta=\sigma$, the
rotation angle is $\pi/2$ making the spin aligned on +$x$ axis that the
precession is maximally suppressed with the oscillation amplitude of the
photocurrent nearly a constant. Nevertheless for $\Delta=0$ and $16$
corresponding to the rotation angle $\phi=0$ and $\pi$, the precession radius
reaches the maximum as well the oscillation amplitude of photocurrent, with a
difference of phase $\pi$.

In summary, by accurately solving the extended Anderson impurity model in the
real time domain with the hierarchical equations of motion (HEOM) approach, we
precisely simulate the whole process of single-hole spin control including
initialization, SU(2) rotation, and read-out. Our theoretical results are in
well agreement with the recent experimental observations, which demonstrates
the feasibility and accuracy for the HEOM approach to describe the hole spin
dynamics. Particularly, the influence of the hybridization strength to
electrodes is fully considered, and a maximal fidelity in the initialization
is predicted.

The support from the NSF of China (No.11374363) and the Research Funds of
Renmin University of China (Grant No. 11XNJ026) is gratefully appreciated.

## References

  * Brunner  _et al._ (2009) D. Brunner, B. D. Gerardot, P. A. Dalgarno, G. Wüst, K. Karrai, N. G. Stoltz, P. M. Petroff, and R. J. Warburton, Science 325, 70 (2009).
  * Press  _et al._ (2008) D. Press, T. D. Ladd, B. Zhang, and Y. Yamamoto, Nature 456, 218 (2008).
  * De Greve  _et al._ (2011) K. De Greve, P. L. McMahon, D. Press, T. D. Ladd, D. Bisping, C. Schneider, M. Kamp, L. Worschech, S. Höfling, A. Forchel, _et al._ , Nat. Phys. 7, 872 (2011).
  * Vamivakas  _et al._ (2010) A. Vamivakas, C.-Y. Lu, C. Matthiesen, Y. Zhao, S. Fält, A. Badolato, and M. Atatüre, Nature 467, 297 (2010).
  * Delteil  _et al._ (2014) A. Delteil, W.-b. Gao, P. Fallahi, J. Miguel-Sanchez, and A. Imamoğlu, Phys. Rev. Lett. 112, 116802 (2014).
  * Vink  _et al._ (2009) I. T. Vink, K. C. Nowack, F. H. Koppens, J. Danon, Y. V. Nazarov, and L. M. Vandersypen, Nat. Phys. 5, 764 (2009).
  * Fischer  _et al._ (2008) J. Fischer, W. Coish, D. Bulaev, and D. Loss, Phys. Rev. B 78, 155329 (2008).
  * Fallahi, Yılmaz, and Imamoğlu (2010) P. Fallahi, S. Yılmaz, and A. Imamoğlu, Phys. Rev. Lett. 105, 257402 (2010).
  * Chekhovich  _et al._ (2011) E. Chekhovich, A. Krysa, M. Skolnick, and A. Tartakovskii, Phys. Rev. Lett. 106, 027402 (2011).
  * Ardelt  _et al._ (2015) P.-L. Ardelt, T. Simmet, K. Müller, C. Dory, K. Fischer, A. Bechtold, A. Kleinkauf, H. Riedl, and J. Finley, Phys. Rev. B 92, 115306 (2015).
  * Gerardot  _et al._ (2008) B. D. Gerardot, D. Brunner, P. A. Dalgarno, P. Öhberg, S. Seidl, M. Kroner, K. Karrai, N. G. Stoltz, P. M. Petroff, and R. J. Warburton, Nature 451, 441 (2008).
  * Godden  _et al._ (2012) T. Godden, J. Quilter, A. Ramsay, Y. Wu, P. Brereton, S. Boyle, I. Luxmoore, J. Puebla-Nunez, A. Fox, and M. Skolnick, Phys. Rev. Lett. 108, 017402 (2012).
  * Godden  _et al._ (2010) T. M. Godden, S. J. Boyle, A. J. Ramsay, A. Fox, and M. Skolnick, arXiv preprint arXiv:1006.5842 (2010).
  * Brash  _et al._ (2015) A. J. Brash, L. M. Martins, F. Liu, J. H. Quilter, A. J. Ramsay, M. S. Skolnick, and A. M. Fox, Phys. Rev. B 92, 121301 (2015).
  * Heiss  _et al._ (2008) D. Heiss, V. Jovanov, M. Bichler, G. Abstreiter, and J. Finley, Phys. Rev. B 77, 235442 (2008).
  * Breuer and Petruccione (2002) H.-P. Breuer and F. Petruccione,  _The theory of open quantum systems_ (Oxford University Press on Demand, 2002).
  * Jin, Zheng, and Yan (2008) J. Jin, X. Zheng, and Y. Yan, J. Chem. Phys. 128, 234703 (2008).
  * Li  _et al._ (2012) Z. Li, N. Tong, X. Zheng, D. Hou, J. Wei, J. Hu, and Y. Yan, Phys. Rev. Lett. 109, 266403 (2012).
  * Ye  _et al._ (2016) L. Ye, X. Wang, D. Hou, R.-X. Xu, X. Zheng, and Y. Yan, WIREs. Comput. Mol. Sci. (2016).
  * Economou and Reinecke (2007) S. E. Economou and T. Reinecke, Phys. Rev. Lett. 99, 217401 (2007).
  * Mar  _et al._ (2013) J. Mar, J. Baumberg, X. Xu, A. Irvine, C. Stanley, and D. Williams, Phys. Rev. B 87, 155315 (2013).
  * Poyatos, Cirac, and Zoller (1997) J. Poyatos, J. I. Cirac, and P. Zoller, Phys. Rev. Lett. 78, 390 (1997).
  * Economou  _et al._ (2006) S. E. Economou, L. Sham, Y. Wu, and D. Steel, Phys. Rev. B 74, 205415 (2006).
  * Kim  _et al._ (2010) E. D. Kim, K. Truex, X. Xu, B. Sun, D. Steel, A. Bracker, D. Gammon, and L. Sham, Phys. Rev. Lett. 104, 167401 (2010).
  * Greilich  _et al._ (2009) A. Greilich, S. E. Economou, S. Spatzek, D. Yakovlev, D. Reuter, A. Wieck, T. Reinecke, and M. Bayer, Nat. Phys. 5, 262 (2009).


# The Derivation of The Probability of Error for BPSK, 16-QAM and 64-QAM in
Rayleigh Fading Channel: A Unified Approach

1st Amr Abdelbari  Artificial Intelligence Engineering Department, AI and
Robotics Institute
Near East University
Nicosia, TRNC, Mersin 10, Turkey
<EMAIL_ADDRESS>2nd Bülent Bilgehan  Dept. of Electrical and
Electronic
Near East University
Nicosia, TRNC, Mersin 10, Turkey
<EMAIL_ADDRESS>
###### Abstract

Understanding the probability of error is paramount in the design and analysis
of digital communication systems, particularly in Rayleigh fading channels
where signal impairments are prevalent. This article presents a unified
approach for deriving the probability of error formulations applicable to
Binary Phase Shift Keying (BPSK), 16-Quadrature Amplitude Modulation (16-QAM),
and 64-QAM in Rayleigh fading channels. This article presents a general
approach to derive the probability of error formulations applicable to Binary
Phase Shift Keying (BPSK), 16-Quadruple Amplitude Modulation (16-QAM) and
64-QAM in Rayleigh fading channels. The derivation process encompasses the
unique characteristics of each modulation scheme and the statistical
properties of Rayleigh fading providing a comprehensive framework to analyze
error performance. By establishing a unified formulation, this approach
simplifies the analysis and facilitates a deeper understanding of error
probability behavior across different modulation schemes. The derived
formulations offer insights into the impact of fading-induced impairments on
system performance, allowing for accurate prediction and optimization of
communication systems in real-world fading environments. The insights and
techniques presented herein serve as a valuable resource for researchers,
engineers, and practitioners engaged in the design and optimization of
communication systems operating in challenging fading environments.

###### Index Terms:

Bit error rate (BER), probability of error, BPSK, 16-QAM, 64-QAM

##  I Introduction

In the realm of digital communication systems, the reliable transmission of
information is paramount, especially in environments characterized by fading
channels. Rayleigh fading, a prevalent phenomenon in wireless communication,
introduces random variations in the received signal amplitude and phase due to
multipath propagation and environmental factors. Understanding and mitigating
the effects of Rayleigh fading are essential for ensuring robust and efficient
communication systems.

Modulation schemes play a crucial role in shaping the performance of
communication systems, offering different trade-offs between spectral
efficiency, complexity, and robustness to channel impairments. Among the
various modulation techniques, BPSK, 16-QAM, and 64-QAM are widely utilized in
practical applications due to their simplicity and efficiency.

However, the integration of these modulation schemes into Rayleigh fading
channels poses significant challenges, primarily due to the unpredictable
nature of fading-induced impairments. Analyzing and quantifying the
performance of BPSK, 16-QAM, and 64-QAM in fading channels requires a thorough
understanding of the probability of error, which measures the likelihood of
incorrect symbol detection at the receiver.

In this paper, we present a unified approach to derive the probability of
error formulations for BPSK, 16-QAM, and 64-QAM in Rayleigh fading channels.
Our approach encompasses the statistical properties of Rayleigh fading and the
modulation characteristics of each scheme, providing a comprehensive framework
for analyzing error performance. By establishing a unified formulation, we aim
to simplify the analysis process and facilitate a deeper understanding of
error behavior across different modulation schemes.

In this paper, we have consistently used specific notation to clarify and
standardize the descriptions of statistical distributions and properties.
i.i.d. stands for independent identically-distributed, indicating that the
elements or samples discussed are both independent from each other and share
the same probability distribution. $\sim$ is used to denote statistically
distributed as. It precedes the description of the statistical distribution
that a variable follows, linking variables to their respective distributions
clearly and succinctly. $\mathcal{CN}(0,\sigma)$ represents a complex normal
vector. This notation specifies that the vector has a zero mean and a variance
denoted by $\sigma$. The complex normal distribution, often essential in
communications theory, particularly in modeling noise or other random
processes, is fundamental for understanding noise characteristics in various
signal processing contexts.

Figure 1: The probability density function PDF for BPSK in AWGN.

##  II The Calculation of BER for BPSK modulation scheme in AWGN channel

The derivation begins by understanding the properties of BPSK modulation
within an Additive White Gaussian Noise (AWGN) channel. BPSK represents binary
digits, 1 and 0, through analog levels $+\sqrt{E_{b}}$ and $-\sqrt{E_{b}}$
respectively. The key characteristics of AWGN are: additive noise is added to
the signal, not multiplied, white noise has a flat spectrum across all
frequencies and finally, Gaussian noise values follow a Gaussian distribution,
described mathematically as:

$p(x)=\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{(x-\mu)^{2}}{2\sigma^{2}}},$
(1)

with $\mu=0$ and $\sigma^{2}=\frac{N_{0}}{2}$. The received signal is either
$y=s_{1}+n$ when bit 1 is transmitted, or $y=s_{0}+n$ when bit 0 is
transmitted. The conditional probability distribution functions (PDFs) for
these scenarios are:

$p(y|s_{0})=\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{(y+\sqrt{E_{b}})^{2}}{2\sigma^{2}}},$
(2)

and

$p(y|s_{1})=\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{(y-\sqrt{E_{b}})^{2}}{2\sigma^{2}}}.$
(3)

Given that both $s_{0}$ and $s_{1}$ are equally probable, the optimal decision
boundary is at 0. Thus, if $y>0$, $s_{1}$ is assumed transmitted; if $y\leq
0$, $s_{0}$ is assumed transmitted.

####  II-1 Probability of error given $s_{1}$ was transmitted

For $s_{1}$ transmission, the probability of error is determined by the
integral of the tail of the Gaussian distribution extending from $-\infty$ to
0:

$p(e|s_{1})=\frac{1}{\sqrt{2\pi\sigma^{2}}}\int_{-\infty}^{0}e^{-\frac{(y-\sqrt{E_{b}})^{2}}{2\sigma^{2}}}dy.$
(4)

Using a substitution where $z^{2}=\dfrac{(y-\sqrt{E_{b}})^{2}}{N_{0}}$ and
$N_{0}=2\sigma^{2}$. Thus, $z=\dfrac{y-\sqrt{E_{b}}}{\sqrt{N_{0}}}$,
$dz=\dfrac{1}{\sqrt{N_{0}}}dy$, and $dy=\sqrt{N_{0}}dz$. The integration
limits will be as follows: $y=0\Rightarrow z=-\sqrt{\dfrac{E_{b}}{N_{0}}}$ and
$y=\infty\Rightarrow z=\infty$. Equation (4) simplifies to:

$\begin{split}p(e|s_{1})=&\dfrac{\sqrt{N_{0}}}{\sqrt{\pi
N_{0}}}\int_{\sqrt{\frac{E_{b}}{N_{0}}}}^{\infty}e^{-z^{2}}dz\\\
=&\dfrac{2}{2\sqrt{\pi}}\int_{\sqrt{\frac{E_{b}}{N_{0}}}}^{\infty}e^{-z^{2}}dz\\\
=&\dfrac{1}{2}erfc(\sqrt{\frac{E_{b}}{N_{0}}}),\end{split}$ (5)

where

$erfc(x)=\dfrac{2}{\sqrt{\pi}}\int_{x}^{\infty}e^{-x^{2}}dx.$ (6)

####  II-2 Probability of error given $s_{0}$ was transmitted

Similarly, the error probability given $s_{0}$ was transmitted is also:

$\begin{split}p(e|s_{0})=&\dfrac{1}{\sqrt{2\pi\sigma^{2}}}\int_{0}^{\infty}e^{\dfrac{-(y+\sqrt{E_{b}})^{2}}{2\sigma^{2}}}dy\\\
=&\dfrac{1}{\sqrt{\pi}}\int_{\sqrt{\frac{E_{b}}{N_{0}}}}^{\infty}e^{-z^{2}}dy\\\
=&\dfrac{1}{2}erfc(\sqrt{\frac{E_{b}}{N_{0}}}).\end{split}$ (7)

####  II-3 Total probability of bit error

Generally, the total probability of bit error can be calculated as follows:

$P_{b}=p(s_{1})p(e|s_{1})+p(s_{0})p(e|s_{0}).$ (8)

Considering equal probabilities for $s_{0}$ and $s_{1}$, i.e.
$p(s_{0})=p(s_{1})=\frac{1}{2}$, the total probability of bit error is:

$P_{b}=\dfrac{1}{2}erfc\big{(}\sqrt{\frac{E_{b}}{N_{0}}}\big{)}.$ (9)

##  III The Calculation of BER for BPSK modulation scheme in Rayleigh Channel

When extending the analysis to a Rayleigh fading channel, the channel response
$h=h_{re}+jh_{im},$ is modeled as a circularly symmetric complex Gaussian
random variable, where the magnitude $|h|$ follows a Rayleigh distribution.
This channel model reflects an environment with numerous reflectors. The real
$h_{re}$ and imaginary $h_{im}$ parts are i.i.d. with mean 0 and variance
$\sigma_{2}$. The probability density of the magnitude $|h|$ is,

$p(h)=\dfrac{h}{\sqrt{\sigma^{2}}}e^{\dfrac{-h^{2}}{2\sigma^{2}}}.$ (10)

###  III-A System model

The received signal in a Rayleigh fading channel is:

$y=hx+n,$ (11)

where $x$ represents the transmitted symbol ($\pm 1$) and $n$ is the AWGN. The
channel is assumed to be flat fading and randomly varying in time, with the
channel $h$ known at the receiver, allowing for straightforward equalization.
The equalization technique employed in this scenario involves a
straightforward division of the received signal $y$ by the apriori known
channel coefficient $h$.the estimated transmitted symbol after equalization is

$\hat{y}=\frac{y}{h}=\frac{hx+n}{h}=x+\frac{n}{h}=x+\tilde{n},$ (12)

where $\tilde{n}=\frac{n}{h}$ denotes the scaled noise, resulting from the
division of the original noise by the channel coefficient.

###  III-B Bit Error Rate

In a Rayleigh fading scenario, the effective bit energy to noise ratio becomes
$\frac{|h|^{2}E_{b}}{N_{0}}$. The conditional BER given the channel $h$ is:

$P_{b|h}=\dfrac{1}{2}erfc\big{(}\sqrt{\frac{|h|^{2}E_{b}}{N_{0}}}\big{)}=\dfrac{1}{2}erfc\big{(}\sqrt{\gamma}\big{)},$
(13)

where $\gamma=\frac{|h|^{2}E_{b}}{N_{0}}$. The unconditional BER incorporates
the distribution of $|h|^{2}$, which affects the overall performance due to
the fading nature of the channel [1].

It is known that if $|h|$ is a Rayleigh distributed random variable, then
$|h|^{2}$ is chi-square distributed with two degrees of freedom. since
$|h|^{2}$ is chi square distributed, $\gamma$ is also chi square distributed.
The probability density function of $\gamma$ is,

$\begin{split}p(\gamma)=&\frac{1}{E_{b}/N_{0}}e^{\frac{-\gamma}{E_{b}/N_{0}}},~{}\gamma\geq
0\\\ =&\frac{1}{\bar{\gamma}}e^{\frac{-\gamma}{\bar{\gamma}}}\\\ \end{split}$
(14)

where $\bar{\gamma}=E_{b}/N_{0}$.

To determine the overall error probability in a Rayleigh fading channel, we
integrate the conditional probability of error over all possible values of the
channel power $\gamma$. The probability of bit error $P_{b}$ in a Rayleigh
fading channel is computed by

$P_{b}=\int_{0}^{\infty}\dfrac{1}{2}erfc\big{(}\sqrt{\gamma}\big{)}p(\gamma)d\gamma.$
(15)

To evaluate the integral, we leverage several mathematical properties and
functions detailed in the Appendix A, and other common mathematical resources.
Specifically, the derivatives of the error function (erf) and its complement
(erfc) are used to facilitate the integration by parts, a technique employed
to simplify the evaluation of integrals involving product terms. The
derivative of the complementary error function with respect to a variable $x$
is given by:

$\dfrac{d}{dx}erfc\big{(}\sqrt{x}\big{)}=-\dfrac{1}{\sqrt{\pi}}e^{-x^{2}}x^{-\frac{1}{2}}$
(16)
$\dfrac{d}{dx}erf\big{(}\sqrt{x}\big{)}=\dfrac{1}{\sqrt{\pi}}e^{-x^{2}}x^{-\frac{1}{2}}$
(17)

These derivatives are crucial for handling the exponential terms that arise
during the integration by parts process. The method of integration by parts is
applied to solve the integral involving $erfc\big{(}\sqrt{\gamma}\big{)}$ and
the exponential decay of $\gamma$ as described by its probability distribution
function (PDF). The integral setup typically follows the formula ( $\int
udv=uv-\int vdu$ ), where choosing appropriate functions for $u$ and $dv$ is
essential to simplify the integration process. For our case, $u$ might be
chosen as $erfc(\sqrt{x})$, $du=-\dfrac{1}{\sqrt{x\pi}}e^{-x}dx$ and $dv$ as
$e^{-\frac{x}{\alpha}}dx$, $v=-\alpha e^{-\frac{x}{\alpha}}$ allowing us to
express the integral in terms that can be more straightforwardly evaluated or
approximated.

$\begin{split}\int erfc\big{(}\sqrt{x}\big{)}e^{-\frac{x}{\alpha}}dx=&-\alpha
erfc(x)e^{-\frac{x}{\alpha}}-\\\
&\int(-\alpha)e^{-\frac{x}{\alpha}}(-\frac{1}{\sqrt{x\pi}})e^{-x}dx\\\
\end{split}$ (18)

The second part can be solved as follows:

$\begin{split}\frac{\alpha}{\sqrt{\pi}}\int
e^{-\frac{x}{\alpha}}e^{-x}x^{-\frac{1}{2}}dx=&\frac{\alpha}{\sqrt{\pi}}\int
e^{-\frac{x(\alpha+1)}{\alpha}}x^{-\frac{1}{2}}dx\\\ \end{split}$ (19)

By taking $u=\frac{x(\alpha+1)}{\alpha}$, $x=\frac{(\alpha)}{\alpha+1}u$, and
$dx=\frac{(\alpha)}{\alpha+1}du$, the previous integral reduces to:

$\begin{split}\frac{\alpha}{\sqrt{\pi}}\int
e^{-\frac{x}{\alpha}}e^{-x}x^{-\frac{1}{2}}dx=&\frac{\alpha}{\sqrt{\pi}}\int
e^{-u}(\frac{\alpha}{\alpha+1})^{-\frac{1}{2}}u^{-\frac{1}{2}}(\frac{\alpha}{\alpha+1})du\\\
=&\frac{\alpha}{\sqrt{\pi}}(\frac{\alpha}{\alpha+1})^{\frac{1}{2}}\int
e^{-u}u^{-\frac{1}{2}}du\\\ \end{split}$ (20)

By taking $t=\sqrt{u}$, $t^{2}=u$, and $du=2tdt$, the previous integral
reduces to:

$\begin{split}\frac{\alpha}{\sqrt{\pi}}(\frac{\alpha}{\alpha+1})^{\frac{1}{2}}\int
e^{-u}u^{-\frac{1}{2}}du=&\frac{\alpha}{\sqrt{\pi}}(\frac{\alpha}{\alpha+1})^{\frac{1}{2}}\int
e^{-t^{2}}t^{-\frac{2}{2}}2tdt\\\
=&\frac{2\alpha}{\sqrt{\pi}}(\frac{\alpha}{\alpha+1})^{\frac{1}{2}}\int
e^{-t^{2}}dt\\\
=&\alpha(\frac{\alpha}{\alpha+1})^{\frac{1}{2}}erf(\sqrt{u})\\\
=&\alpha(\frac{\alpha}{\alpha+1})^{\frac{1}{2}}erf(\sqrt{\frac{(\alpha+1)}{\alpha}x})\\\
\end{split}$ (21)

By replacing the last part with the second part in (18):

$\begin{split}\int erfc\big{(}\sqrt{x}\big{)}e^{-\frac{x}{\alpha}}dx=&-\alpha
erfc(x)e^{-\frac{x}{\alpha}}-\\\
&\alpha(\frac{\alpha}{\alpha+1})^{\frac{1}{2}}erf(\sqrt{\frac{(\alpha+1)}{\alpha}x})\\\
\end{split}$ (22)

With the preceding analysis, deriving the BER for BPSK modulation in a
Rayleigh fading channel becomes straightforward. In (22), let
$\alpha=\bar{\gamma}$ and $x=\gamma$:

$\begin{split}P_{b}=&\int_{0}^{\infty}\dfrac{1}{2\bar{\gamma}}erfc\big{(}\sqrt{\gamma}\big{)}e^{\frac{-\gamma}{\bar{\gamma}}}d\gamma\\\
=&\dfrac{1}{2\bar{\gamma}}\big{[}-\bar{\gamma}erfc(\sqrt{\gamma})e^{-\frac{\gamma}{\bar{\gamma}}}-\\\
&\bar{\gamma}\sqrt{\frac{\bar{\gamma}}{\bar{\gamma}+1}}erf(\sqrt{\frac{\bar{\gamma}+1}{\bar{\gamma}}\gamma})\big{]}_{0}^{\infty}\\\
=&\dfrac{1}{2\bar{\gamma}}\big{[}\bar{\gamma}erfc(\sqrt{\gamma})e^{-\frac{\gamma}{\bar{\gamma}}}+\\\
&\bar{\gamma}\sqrt{\frac{\bar{\gamma}}{\bar{\gamma}+1}}erf(\sqrt{\frac{\bar{\gamma}+1}{\bar{\gamma}}\gamma})\big{]}_{\infty}^{0}\\\
\end{split}$ (23)

Taking $erfc(0)=1$, $erfc(\infty)=0$, $erf(0)=0$ and $erf(\infty)=1$:

$\begin{split}P_{b}=&\dfrac{1}{2\bar{\gamma}}\big{[}(\bar{\gamma}+0)-(0+\bar{\gamma}\sqrt{\frac{\bar{\gamma}}{\bar{\gamma}+1}})\big{]}\\\
=&\dfrac{\bar{\gamma}}{2\bar{\gamma}}\big{[}(1-\sqrt{\frac{\bar{\gamma}}{\bar{\gamma}+1}})\big{]}\\\
=&\dfrac{1}{2}\big{[}(1-\sqrt{\frac{E_{b}/N_{0}}{E_{b}/N_{0}+1}})\big{]}\\\
\end{split}$ (24)

Figure 2: Figure: Constellation plot for 4 PAM modulation.

##  IV The Calculation of BER for 4PAM Modulation Scheme in AWGN Channel

Consider the symbols used in a 4-PAM system, which are defined as
$\alpha_{4QAM}={\pm 1,\pm 3}$. The average energy of the PAM constellation,
assuming each symbol is equally likely to occur, can be calculated as follows:

$E_{4PAM}=E\big{[}|\alpha_{4PAM}|^{2}\big{]}=\frac{1}{4}\sum[1^{2}+3^{2}]=5$
(25)

This expression considers the squared amplitude of each symbol, reflecting the
energy per symbol in the constellation [2]. A normalized plot of the 4-PAM
signal’s constellation can be visualized as depicted in Fig. 2. It is assumed
that the additive noise, $n$, adhering to the Gaussian probability
distribution function, impacts the received signals, as described in equation
(1). The received signal for a transmitted symbol $s_{i}$ can be modeled by:

$y=s_{i}+n,\quad i=1,2,3,4$ (26)

Focusing on the scenario where $s_{3}$ is transmitted, the conditional PDF for
the received signal $y$ is:

$p(y|s_{3})=\frac{1}{\sqrt{\pi
N_{0}}}e^{-\frac{(y-3\sqrt{\frac{E_{s}}{5}})^{2}}{N_{0}}}$ (27)

This scenario assumes the detection threshold between the signals for +1 and
+3 is set at $+2\sqrt{\frac{E_{s}}{5}}$. Consequently, if the received signal
$y$ exceeds this threshold, the receiver interprets that $s_{3}$ was
transmitted.

Using this threshold, the probability of error when $s_{3}$ is transmitted -
represented by the area under the curve in the specified region (the blue
area) - can be calculated to ensure accurate detection and error assessment in
noisy environments by

$p(e|s_{3})=\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{-\infty}^{+2\sqrt{\frac{E_{s}}{5}}}e^{-\dfrac{(y-3\sqrt{\frac{E_{s}}{5}})^{2}}{N_{0}}}dy$
(28)

By taking $z^{2}=\dfrac{(y-3\sqrt{\frac{E_{s}}{5}})^{2}}{N_{0}}$. Thus,
$z=\dfrac{y-3\sqrt{\frac{E_{s}}{5}}}{\sqrt{N_{0}}}$,
$dz=\dfrac{1}{\sqrt{N_{0}}}dy$, and $dy=\sqrt{N_{0}}dz$. The integration
limits will be as follows: $y=+2\sqrt{\frac{E_{s}}{5}}\Rightarrow
z=-\sqrt{\dfrac{E_{s}}{5N_{0}}}$ and $y=-\infty\Rightarrow z=-\infty$.
Therefore, (28)will be:

$\begin{split}p(e|s_{3})=&\dfrac{\sqrt{N_{0}}}{\sqrt{\pi
N_{0}}}\int_{\sqrt{\frac{E_{s}}{5N_{0}}}}^{\infty}e^{-z^{2}}dz\\\
=&\dfrac{2}{2\sqrt{\pi}}\int_{\sqrt{\frac{E_{s}}{5N_{0}}}}^{\infty}e^{-z^{2}}dz\\\
=&\dfrac{1}{2}erfc(\sqrt{\frac{E_{s}}{5N_{0}}}).\end{split}$ (29)

[1] Given the symmetrical nature of the 4-PAM constellation, particularly
regarding the symbols +3 and -3, it is intuitively reasonable to conclude that
the probability of error when $s_{3}$ (representing +3) is transmitted mirrors
that when -3 is transmitted. This symmetry in the constellation layout implies
that the error characteristics for these signals are identical due to their
equidistant spacing from their respective decision boundaries in the presence
of Gaussian noise.

Focusing now on the symbol $s_{2}$, which typically represents one of the
intermediate values in a PAM constellation, the probability of error when
$s_{2}$ is transmitted (the area in green and red region) can be calculated by

$\begin{split}p(e|s_{2})=&\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{-\infty}^{0}e^{-\dfrac{(y-\sqrt{\frac{E_{s}}{5}})^{2}}{N_{0}}}dy+\\\
&\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{+2\sqrt{\frac{E_{s}}{5}}}^{\infty}e^{-\dfrac{(y-\sqrt{\frac{E_{s}}{5}})^{2}}{N_{0}}}dy\end{split}$
(30)

For the first two part of (30), take
$z^{2}=\dfrac{(y-\sqrt{\frac{E_{s}}{5}})^{2}}{N_{0}}$. Thus,
$z=\dfrac{y-\sqrt{\frac{E_{s}}{5}}}{\sqrt{N_{0}}}$,
$dz=\dfrac{1}{\sqrt{N_{0}}}dy$, and $dy=\sqrt{N_{0}}dz$. The integration
limits of the first integral will be $y=0\Rightarrow
z=-\sqrt{\dfrac{E_{s}}{5N_{0}}}$ and $y=-\infty\Rightarrow z=-\infty$. The
integration limits of the second integral will be
$y=2\sqrt{\dfrac{E_{s}}{5N_{0}}}\Rightarrow z=\sqrt{\dfrac{E_{s}}{5N_{0}}}$
and $y=\infty\Rightarrow z=\infty$. Therefore, (30) will be:

$\begin{split}p(e|s_{2})=&\dfrac{1}{\sqrt{\pi}}\int_{-\infty}^{-\sqrt{\frac{E_{s}}{5N_{0}}}}e^{-z^{2}}dz+\\\
&\dfrac{1}{\sqrt{\pi}}\int_{\sqrt{\frac{E_{s}}{5N_{0}}}}^{\infty}e^{-z^{2}}dz\\\
=&\dfrac{2}{\sqrt{\pi}}\int_{\sqrt{\frac{E_{s}}{5N_{0}}}}^{\infty}e^{-z^{2}}dz\\\
=&erfc(\sqrt{\frac{E_{s}}{5N_{0}}}).\end{split}$ (31)

Considering the symmetrical arrangement of the +1 and -1 symbols in the 4-PAM
constellation, it logically follows that the probability of error when
transmitting $s_{1}$ (representing +1) is identical to that when transmitting
$s_{2}$ ( representing -1). Since each symbol in the constellation is assumed
to be equally likely, the total probability of symbol error across the
constellation can be determined by

$P_{s}=p(s_{0})p(e|s_{0})+p(s_{1})p(e|s_{1})+p(s_{2})p(e|s_{2})+p(s_{3})p(e|s_{3})$
(32)

Given that $s_{0}$, $s_{1}$, $s_{2}$ and $s_{3}$ are equally probable i.e.
$p(s_{0})=p(s_{1})=p(s_{2})=p(s_{3})=\frac{1}{4}$, the symbol error
probability is,

$P_{s}=\dfrac{1}{4}\big{[}erfc\big{(}\sqrt{\frac{E_{s}}{5N_{0}}}\big{)}+2erfc\big{(}\sqrt{\frac{E_{s}}{5N_{0}}}\big{)}\big{]}=\dfrac{3}{4}erfc\big{(}\sqrt{\frac{E_{s}}{5N_{0}}}\big{)}$
(33)

Since $E_{s}/N_{0}=2E_{b}/N_{0}$ and BER = $\frac{1}{log_{2}(M)}$ SER, the bit
error probability is [3]:

$P_{b}=\dfrac{3}{8}erfc\big{(}\sqrt{\frac{2E_{b}}{5N_{0}}}\big{)}$ (34)

To be expressed in terms of distance between each successive symbols
$d=2\sqrt{\frac{E_{s}}{N_{0}}}$.

$P_{b}=\dfrac{3}{8}erfc\big{(}\frac{d}{2\sqrt{2}\sigma}\big{)}$ (35)

or generally for any M-PAM [4]:

$P_{b}=\dfrac{M-1}{Mlog_{2}(M)}erfc\big{(}\frac{d}{2\sqrt{2}\sigma}\big{)}$
(36)

Figure 3: Constellation plot for QPSK (4-QAM) constellation.

##  V The Calculation of BER for 4-QAM Modulation Scheme in AWGN Channel

Consider that the alphabets used for a QPSK (4-QAM) is $\alpha_{4QAM}={\pm
1\pm 1j}$. The scaling factor of $\sqrt{\frac{E_{s}}{2}}$ is utilized to
normalize the average energy of the transmitted symbols to 1. This
normalization is based on the assumption that all points in the constellation
are equally probable. Such scaling ensures that the total power of the
constellation adheres to the system’s energy constraints, thereby optimizing
the transmission power across different communication scenarios. It is assumed
that the additive noise, denoted by $n$, adheres to the Gaussian probability
distribution function, as specified in (1).

For the symbol $s_{2}$, the conditional PDF when $s_{2}$ is transmitted is
given by:

$p(y|s_{2})=\frac{1}{\sqrt{N_{0}\pi}}e^{-\frac{(y-\sqrt{\frac{E_{s}}{2}})^{2}}{N_{0}}}.$
(37)

This equation indicates how the received symbol $y$ is distributed around the
transmitted symbol $s_{2}$ under Gaussian noise. As illustrated in Fig. 3,
correct decoding of symbol $s_{2}$ is contingent upon $y$ landing within a
specific hashed region. The probability that the symbol $s_{2}$ is correctly
decoded, considering both the real and imaginary components of $y$, is:

$p(c|s_{2})=p(\Re_{y}>0|s_{2})\cdot p(\Im_{y}>0|s_{2})$ (38)

This probability effectively represents the likelihood that the real component
of $y$ is greater than 0 given that $s_{2}$ was transmitted, and is determined
by excluding the area outside the designated red region. The probability that
the real component of the received signal $y$ is greater than 0, given
$s_{2}$was transmitted, is critical for successful demodulation. This
probability, representing the area outside the red region, is given by:

$\begin{split}p(\Re_{y}>0|s_{2})=&1-\dfrac{1}{\sqrt{N_{0}\pi}}\int_{-\infty}^{0}e^{-\dfrac{(y-\sqrt{\frac{E_{s}}{2}})^{2}}{N_{0}}}dy\\\
=&1-erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}\end{split}$ (39)

This result is derived by calculating the area under the Gaussian curve from
negative infinity to zero, subtracted from 1, which indicates the cumulative
probability from zero to positive infinity. Similarly, the probability that
the imaginary component of $y$ is greater than 0, given $s_{2}$ was
transmitted (i.e., the area outside the blue region), is crucial for
determining the accurate decoding of the signal. It is calculated as:

$\begin{split}p(\Im_{y}>0|s_{2})=&1-\dfrac{1}{\sqrt{N_{0}\pi}}\int_{-\infty}^{0}e^{-\dfrac{(y-\sqrt{\frac{E_{s}}{2}})^{2}}{N_{0}}}dy\\\
=&1-erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}.\end{split}$ (40)

The overall probability that $s_{2}$ is correctly decoded, which considers
both the real and imaginary components, is then given by:

$\begin{split}p(c|s_{2})=&\big{[}1-erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}\big{]}^{2}\\\
=&1-\dfrac{2}{2}erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}+\dfrac{1}{4}erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}^{2}\\\
=&1-erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}+\dfrac{1}{4}erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}^{2}\\\
\end{split}$ (41)

The probability that the symbol $s_{2}$ is in error, i.e., at least one
component is incorrectly decoded, is thus:

$\begin{split}P_{QPSK}=&1-p(c|s_{2})\\\
=&1-\big{[}1-erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}+\dfrac{1}{4}erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}^{2}\big{]}\\\
=&erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}-\dfrac{1}{4}erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}^{2}\\\
\end{split}$ (42)

For high values of the signal-to-noise ratio $\frac{E_{s}}{2N_{0}}$, where the
system approaches an ideal state with minimal noise, the second term becomes
negligible, allowing the error probability to be approximated simply as:

$P_{QPSK}\approx erfc\big{(}\sqrt{\frac{E_{s}}{2N_{0}}}\big{)}$ (43)

This approximation highlights that as the SNR increases, the system’s
performance significantly improves, reducing the error rate and enhancing
reliable communication [2].

Figure 4: Constellation plot for 64-QAM modulation (without the scaling factor
of $\frac{1}{\sqrt{42}}$).

##  VI The Calculation of BER for M-QAM Modulation Scheme in AWGN Channel

The M-QAM constellation is defined based on $M=2^{b}$, where $b$ represents
the number of bits per symbol. For practical and analytical convenience, $b$
is often chosen to be even due to the following advantages. First, half of the
bits are mapped onto the real axis and the other half onto the imaginary axis.
This division results in two independent $b/2$-level PAM signals, simplifying
the design of the mapping scheme. Second, decoding processes can be
independently applied to the real and imaginary components of the signal,
which streamlines the receiver architecture. It is important to note, however,
that this square constellation layout is not always the most efficient in
terms of signal-to-noise ratio performance.

For an even $b$ within an M-QAM system, the constellation points are typically
given by:

$\alpha_{MQAM}={\pm(2m-1)\pm(2m-1)j},~{}m\in{1,2,\dots,\frac{\sqrt{M}}{2}}.$
(44)

For instance, a 64-QAM constellation, represented by ($M=64$) constellation,
$m\in{1,2,3,4}$.

$\begin{split}\alpha_{64QAM}=\begin{bmatrix}(\pm 7\pm 7j)&(\pm 7\pm 5j)&(\pm
7\pm 3j)&(\pm 7\pm 1j)\\\ (\pm 5\pm 7j)&(\pm 5\pm 5j)&(\pm 5\pm 3j)&(\pm 5\pm
1j)\\\ (\pm 3\pm 7j)&(\pm 3\pm 5j)&(\pm 3\pm 3j)&(\pm 3\pm 1j)\\\ (\pm 1\pm
7j)&(\pm 1\pm 5j)&(\pm 1\pm 3j)&(\pm 1\pm 1j)\\\ \end{bmatrix}.\end{split}$

(45)

To calculate the average energy across the constellation, calculate the
cumulative energy of individual symbols.

$E_{4\alpha}=\sum_{m=1}^{\frac{\sqrt{M}}{2}}\big{[}(2m-1)+j(2m-1)\big{]}^{2}=\frac{\sqrt{M}}{2}(M-1).$
(46)

Each symbol appears $2\sqrt{M}$ times within the constellation. The average
energy is then:

$E_{MQAM}=\frac{2\sqrt{M}}{M}E_{\alpha}=\frac{2\sqrt{M}}{M}\frac{\sqrt{M}}{3}(M-1)=\frac{2}{3}(M-1)$
(47)

For example, for 64-QAM, the average energy is

$E_{64-QAM}=\frac{2}{3}(64-1)=42,$ (48)

and for 16-QAM,

$E_{16-QAM}=\frac{2}{3}(16-1)=10.$ (49)

For specific constellations like 64-QAM and 16-QAM, the respective average
energies $\frac{1}{\sqrt{10}}$, $\frac{1}{\sqrt{42}}$ can be calculated and
used to normalize the transmission power [2].

To analyze the symbol error rate, especially within the context of a 64-QAM
system and then extended to general M-QAM formats. Fig. 4 illustrates the
differentiation of constellation points. The corner points shaped as square
(in red) always four points $N_{corner}=4$. The inside points shaped as
diamond (in magenta) computed as $N_{corner}=(\sqrt{M}-2)^{2}$. Last points
are edge points which are neither corners nor center shaped as stars (in
blue). The number of these points is calculated as $N_{side}=4(\sqrt{M}-2)$.
For 64-QAM it is 24 points.

Assuming the received symbol $y$ is affected by $n$ the AWGN in (1), the
received value is:

$y=k\sqrt{E_{s}}+n.$ (50)

Here, $k=\sqrt{\dfrac{1}{\frac{2}{3}(M-1)}}$ acts as the normalizing factor
based on the constellation’s average energy, and $E_{s}$ is the symbol energy.

####  VI-1 For constellation points located inside the M-QAM grid

, particularly for those not on the edges or corners, the conditional PDF for
the received signal $y$ , given that the symbol
${+k\sqrt{E_{s}},+k\sqrt{E_{s}}}$ is transmitted, is modeled as:

$p(y|inside)=\dfrac{1}{\sqrt{N_{0}\pi}}e^{-\dfrac{(y-k\sqrt{E_{s}})^{2}}{N_{0}}}$
(51)

Correct decoding of the symbol ($I=+1,Q=+1$)is contingent upon the real
$y~{}(\Re_{y})$and imaginary $y~{}(\Im_{y})$ components of the received signal
$y$ falling within specific bounds:

$p(c|inside)=p(\Re_{y}>0,\Re_{y}\leq 2k\sqrt{E_{s}}|+1)p(\Im_{y}>0,\Im_{y}\leq
2k\sqrt{E_{s}}|+1)$

(52)

To compute this probability, calculate the probability that the real component
$y~{}(\Re_{y})$ exceeds 0 and does not surpass $2k\sqrt{E_{s}}$, and similarly
for the imaginary component $y~{}(\Im_{y})$ . Then integrate the conditional
PDF from 0 to $2k\sqrt{E_{s}}$ for both the real and imaginary components.

For instance, to determine the likelihood that the real component
$y~{}(\Re_{y})$ of the received signal falls within the specified range of 0
to 2, the calculation involves the integration of the PDF across two distinct
regions: 1) Probability of exceeding the upper limit: calculate the
probability that the real component exceeds the upper limit of 2, extending to
infinity $\infty$. 2) Probability of falling below the lower limit: calculate
the probability that the real component is less than the lower boundary of 0,
extending to negative infinity $-\infty$. Since the total probability for any
variable within its complete range is 1, the probability that the real
component falls within the range from 0 to 2 is obtained by subtracting the
sum of the probabilities calculated in steps 1 and 2 from 1. This approach
uses the complement rule to derive the probability for the desired range as
follows:

$p(\Re_{y}\text{ within }0\text{ to }2)=1-[p(\Re_{y}>2)+p(\Re_{y}<0)]$ (53)

This method effectively segments the total probability distribution into parts
that are outside the interval of interest and subtracts their cumulative
impact from unity to find the probability of the variable lying within the
specified range.

In general, using the CDF and the error function (erfc), the probability for
each component can be derived from:

$\begin{split}p(\Re_{y}>0,\Re_{y}\leq 2k\sqrt{E_{s}}|+1)=&1-\\\
&\big{[}\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{-\infty}^{0}e^{-\dfrac{(y-k\sqrt{E_{s}})^{2}}{N_{0}}}dy+\\\
&\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{+2k\sqrt{E_{s}}}^{\infty}e^{-\dfrac{(y-k\sqrt{E_{s}})^{2}}{N_{0}}}dy\big{]}\\\
=&1-erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\end{split}$

(54)

Similarly,

$p(\Im_{y}>0,\Im_{y}\leq
2k\sqrt{E_{s}}|+1)=1-erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}$ (55)

The combined probability that the symbol is decoded correctly (both real and
imaginary parts fall within the designated range) is:

$p(c|inside)=\big{[}1-erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}\big{[}1-erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}.$
(56)

And consequently, the error probability is:

$\begin{split}p(e|inside)=&1-\big{[}1-erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}^{2}\\\
=&2erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}-erfc^{2}\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}.\\\
\end{split}$ (57)

####  VI-2  For corner symbols in an M-QAM constellation, such as $I=+7,Q=+7$

, the PDF for the received signal $y$, given that this specific symbol
${+7k\sqrt{E_{s}},+7k\sqrt{E_{s}}}$ was transmitted, is defined as:

$p(y|corner)=\dfrac{1}{\sqrt{N_{0}\pi}}e^{-\dfrac{(y-7k\sqrt{E_{s}})^{2}}{N_{0}}}$
(58)

This formula calculates the likelihood of $y$ based on its deviation from the
expected position of the corner symbol, adjusted by the noise power $N_{0}$.
To determine if the corner symbol is correctly decoded, both the real
$y~{}(\Re_{y})$ and imaginary $y~{}(\Im_{y})$ parts of the received symbol
must exceed a threshold of 6. This ensures the symbol is decoded as
$I=+7,Q=+7$ and not mistaken for a neighboring symbol. The probability of
correct decoding is expressed as:

$p(c|corner)=p(\Re_{y}>6,\Re_{y}\leq\infty|+7)p(\Im_{y}>6,\Im_{y}\leq\infty|+7)$

(59)

The thresholds are determined by integrating the tail of the conditional PDF
from 6 to $\infty$. This provides the probability that $y$ exceeds the
threshold, indicating correct symbol detection:

$\begin{split}p(\Re_{y}>6,\Re_{y}\leq\infty|+7)=&\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{6k\sqrt{E_{s}}}^{\infty}e^{-\dfrac{(y-7k\sqrt{E_{s}})^{2}}{N_{0}}}dy\\\
=&1-\dfrac{1}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}.\end{split}$

(60)

Similarly, for the imaginary part:

$p(\Im_{y}>6,\Im_{y}\leq\infty|+7)=1-\dfrac{1}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}.$
(61)

By multiplying the probabilities for the real and imaginary components, the
overall probability of correctly decoding the corner symbol is:

$p(c|corner)=\big{[}1-\dfrac{1}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}\big{[}1-\dfrac{1}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}.$
(62)

The likelihood that the symbol is decoded erroneously (i.e., at least one
component is incorrectly decoded) is given by the complement of the correct
decoding probability:

$\begin{split}p(e|corner)=&1-\big{[}1-\dfrac{1}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}^{2}\\\
=&erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}-\dfrac{1}{4}erfc^{2}\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}.\\\
\end{split}$ (63)

####  VI-3 For symbols that are neither at the corners nor strictly inside the
grid, such as $I=+7,Q=+1$

, a detailed error analysis involves assessing the real and imaginary
components of the received signal $y$. Correct decoding for these symbols is
contingent upon the components falling within specific ranges. The probability
of correct demodulation is,

$p(c|side)=p(\Re_{y}>6,\Re_{y}\leq\infty|+7)p(\Im_{y}>0,\Im_{y}\leq
2k\sqrt{E_{s}}|+1).$

(64)

Correct decoding of this symbol requires the real part of $y~{}(\Re_{y})$ must
be greater than 6 and can extend to $\infty$, and the imaginary part of
$y~{}(\Im_{y})$ must lie between 0 and 2.

The probability that $y~{}(\Re_{y})$ falls within the correct range is
calculated by integrating the tail of the PDF from 6 to $\infty$:

$\begin{split}p(\Re_{y}>6,\Re_{y}\leq\infty|+7)=&\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{6k\sqrt{E_{s}}}^{\infty}e^{-\dfrac{(y-7k\sqrt{E_{s}})^{2}}{N_{0}}}dy\\\
=&1-\dfrac{1}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\end{split}$

(65)

For the imaginary component, the integration considers two segments: 1) From
$\infty$ to 0 to find the lower boundary. 2) From $2k\sqrt{E_{s}}$ to $\infty$
for the upper boundary. The probability that $y~{}(\Im_{y})$ falls within 0 to
2 is:

$\begin{split}p(\Im_{y}>0,\Im_{y}\leq 2k\sqrt{E_{s}}|+1)=&1-\\\
&\big{[}\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{-\infty}^{0}e^{-\dfrac{(y-k\sqrt{E_{s}})^{2}}{N_{0}}}dy+\\\
&\dfrac{1}{\sqrt{\pi
N_{0}}}\int_{+2k\sqrt{E_{s}}}^{\infty}e^{-\dfrac{(y-k\sqrt{E_{s}})^{2}}{N_{0}}}dy\big{]}\\\
=&1-erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}.\end{split}$

(66)

Combining these probabilities, the overall likelihood that the symbol
$I=+7,Q=+1$ is decoded correctly is:

$\begin{split}p(c|side)=&\big{[}1-erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}\big{[}1-\dfrac{1}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\big{]}\\\
=&1-\dfrac{3}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}+\dfrac{1}{2}erfc^{2}\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\\\
\end{split}$ (67)

Finally, the probability that the symbol is decoded incorrectly, which occurs
if either the real or imaginary part falls outside their respective correct
ranges, is:

$\begin{split}p(e|side)=&1-p(c|side)\\\
=&\dfrac{3}{2}erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}+\dfrac{1}{2}erfc^{2}\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}\\\
\end{split}$ (68)

###  VI-A Symbol Error Probabilities Across Constellation Points

Given the previously calculated error probabilities for different types of
constellation points within an M-QAM system, the total SER is computed by
averaging the error probabilities across all constellation point categories.
This approach considers the varying likelihood of error occurrence depending
on the location of each symbol within the constellation grid. The joint symbol
error rate is calculated using the formula:

$\begin{split}P(e|MQAM)=&\dfrac{1}{M}\big{[}N_{inside}p(e|inside)+\\\
&N_{corner}p(e|corner)+N_{side}p(e|side)\big{]}.\end{split}$ (69)

This equation effectively weights the error probabilities by the number of
symbols in each category relative to the total number of symbols $M$.
Integrating the individual probabilities:

$\begin{split}P(e|MQAM)=&2(1-\dfrac{1}{\sqrt{M}})erfc\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}-\\\
&(1-\dfrac{2}{\sqrt{M}}+\dfrac{1}{M})erfc^{2}\big{(}k\sqrt{\frac{E_{s}}{N_{0}}}\big{)}.\\\
\end{split}$ (70)

Utilizing the $Q-$function to express the error probability provides a more
intuitive understanding of the system’s performance under Gaussian noise [5]:

$\begin{split}P(e|MQAM)=&4(\dfrac{\sqrt{M}-1}{\sqrt{M}})Q\big{(}\sqrt{\frac{3E_{s}}{N_{0}(M-1)}}\big{)}-\\\
&4(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}Q^{2}\big{(}\sqrt{\frac{3E_{s}}{N_{0}(M-1)}}\big{)}\\\
\end{split}$ (71)

This representation is particularly useful for communications engineers
designing systems to meet specific performance criteria under normal operating
conditions. The BER is derived from the SER considering the modulation
efficiency. $E_{s}/N_{0}=qE_{b}/N_{0}$, $q=log_{2}(M)$, the BER =
$\frac{1}{log_{2}(M)}$SER. For a 16-QAM system example, SER is calculated as
[6]:

$\begin{split}P_{s}(e|16QAM)=&\dfrac{3}{2}erfc\big{(}\sqrt{\frac{E_{s}}{10N_{0}}}\big{)}-\dfrac{9}{16}erfc^{2}\big{(}\sqrt{\frac{E_{s}}{10N_{0}}}\big{)}\\\
\approx&\dfrac{3}{2}erfc\big{(}\sqrt{\frac{E_{s}}{10N_{0}}}\big{)}\end{split}$
(72)

The corresponding BER for Gray coded 16-QAM in an AWGN environment where $q=4$
is:

$P_{b}(e|16QAM)=\dfrac{3}{8}erfc\big{(}\sqrt{\frac{2E_{b}}{5N_{0}}}\big{)}$
(73)

##  VII The Calculation of BER for M-QAM Modulation Scheme in Rayleigh Fading
Channel

To get the probability of error over Rayleigh fading, we should average the
BER in AWGN over the fading channel as follow [7]:

$P_{b}=\int_{0}^{\infty}P_{e}(\gamma)p_{\gamma}(\gamma)d\gamma,$ (74)

where the probability density function of $\gamma$ is,

$\begin{split}p_{\gamma}(\gamma)=&\frac{1}{E_{b}/N_{0}}e^{\frac{-\gamma}{E_{b}/N_{0}}},~{}\gamma\geq
0\\\ =&\frac{1}{\bar{\gamma}}e^{\frac{-\gamma}{\bar{\gamma}}},\\\ \end{split}$
(75)

and $\bar{\gamma}=E_{b}/N_{0}$. WE will use the following Craig forms for the
Q-function.

$Q(x)=\dfrac{1}{\pi}\int_{0}^{\frac{\pi}{2}}e^{-\dfrac{x^{2}}{2\sin^{2}(\theta)}}d\theta,$
(76)

and

$Q^{2}(x)=\dfrac{1}{\pi}\int_{0}^{\frac{\pi}{4}}e^{-\dfrac{x^{2}}{2\sin^{2}(\theta)}}d\theta.$
(77)

By replacing $x{2}=\dfrac{3\gamma}{M-1}$ in (76) and (77). Then by
substituting them in (71) and (74), the probability of error in Rayleigh
fading will be

$\begin{split}P_{e-MQAM}=&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{0}^{\frac{\pi}{2}}\int_{0}^{\infty}e^{-\frac{3}{M-1}\frac{\gamma}{2\sin^{2}(\theta)}}d\theta
p_{\gamma}(\gamma)d\gamma-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\int_{0}^{\frac{\pi}{4}}\int_{0}^{\infty}e^{-\frac{3}{M-1}\frac{\gamma}{2\sin^{2}(\theta)}}d\theta
p_{\gamma}(\gamma)d\gamma.\\\ \end{split}$ (78)

To solve this we will use the Laplace transformation which implies:

$\begin{split}F_{\gamma}(S)=&\int_{0}^{\infty}e^{-S\gamma}p_{\gamma}(\gamma)d\gamma\\\
=&\int_{0}^{\infty}e^{-S\gamma}\frac{1}{\bar{\gamma}}e^{\frac{-\gamma}{\bar{\gamma}}}d\gamma\\\
=&\dfrac{1}{1+S\bar{\gamma}}.\end{split}$ (79)

Let $S=(\frac{3}{M-1})\frac{1}{2\sin^{2}(\theta)}$.

$\begin{split}P_{e-MQAM}=&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{0}^{\frac{\pi}{2}}F_{\gamma}\big{(}(\frac{3}{M-1})\frac{1}{2\sin^{2}(\theta)}\big{)}d\theta-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\int_{0}^{\frac{\pi}{4}}F_{\gamma}\big{(}(\frac{3}{M-1})\frac{1}{2\sin^{2}(\theta)}\big{)}d\theta\\\
=&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{0}^{\frac{\pi}{2}}\big{(}\dfrac{1}{1+(\frac{3}{M-1})\frac{\bar{\gamma}}{2\sin^{2}(\theta)}}\big{)}d\theta-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\int_{0}^{\frac{\pi}{4}}\big{(}\dfrac{1}{1+(\frac{3}{M-1})\frac{\bar{\gamma}}{2\sin^{2}(\theta)}}\big{)}d\theta\\\
=&P_{e-MQAM}(1)-P_{e-MQAM}(2).\end{split}$ (80)

Let’s solve the first part of the previous equation:

$\begin{split}P_{e-MQAM}(1)=&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{0}^{\frac{\pi}{2}}\big{(}1+(\frac{3}{M-1})\frac{\bar{\gamma}}{2\sin^{2}(\theta)}\big{)}^{-1}d\theta\\\
=&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{0}^{\frac{\pi}{2}}\big{(}1-1+\dfrac{2\sin^{2}(\theta)}{2\sin^{2}(\theta)+(\frac{3}{M-1})\bar{\gamma}}\big{)}d\theta\\\
=&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{0}^{\frac{\pi}{2}}\big{(}1-\dfrac{2\sin^{2}(\theta)+(\frac{3}{M-1})\bar{\gamma}}{2\sin^{2}(\theta)+(\frac{3}{M-1})\bar{\gamma}}+\\\
&\dfrac{2\sin^{2}(\theta)}{2\sin^{2}(\theta)+(\frac{3}{M-1})\bar{\gamma}}\big{)}d\theta\\\
=&2(\dfrac{\sqrt{M}-1}{\sqrt{M}})-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{0}^{\frac{\pi}{2}}\big{(}\dfrac{(\frac{3}{M-1})\bar{\gamma}}{2\sin^{2}(\theta)+(\frac{3}{M-1})\bar{\gamma}}\big{)}d\theta.\\\
\end{split}$ (81)

Now, let $\tan(\theta)=t$,
$\dfrac{dt}{d\theta}=\sec^{2}(\theta)=1+\tan^{2}(\theta)=1+t^{2}$,
$d\theta=\dfrac{1}{1+t^{2}}dt$, and $\sin(\theta)=\dfrac{t}{\sqrt{1+t^{2}}}$.

$\begin{split}P_{e-MQAM}(1)=&2(\dfrac{\sqrt{M}-1}{\sqrt{M}})-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})\int_{\tan(0)}^{\tan(\frac{\pi}{2})}\big{(}\dfrac{(\frac{3}{M-1})\bar{\gamma}}{2\dfrac{t^{2}}{1+t^{2}}+(\frac{3}{M-1})\bar{\gamma}}\big{)}\\\
&*\dfrac{1}{1+t^{2}}dt\\\ =&2(\dfrac{\sqrt{M}-1}{\sqrt{M}})-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})(\frac{3}{M-1})\bar{\gamma}\\\
&*\int_{\tan(0)}^{\tan(\frac{\pi}{2})}\big{(}\dfrac{1+t^{2}}{2t^{2}+(\frac{3}{M-1})\bar{\gamma}(1+t^{2})}\big{)}\\\
&*\dfrac{1}{1+t^{2}}dt\\\ =&2(\dfrac{\sqrt{M}-1}{\sqrt{M}})-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})(\frac{3}{M-1})\bar{\gamma}\\\
&*\int_{\tan(0)}^{\tan(\frac{\pi}{2})}\big{(}\dfrac{1}{(\frac{3}{M-1})\bar{\gamma}+(2+(\frac{3}{M-1})\bar{\gamma})t^{2}}\big{)}dt.\\\
\end{split}$ (82)

To solve the integral, let $a=\frac{3}{M-1})\bar{\gamma}$ and
$b=2+\frac{3}{M-1})\bar{\gamma}$ [6]. Then the integral will be:

$\begin{split}\int_{\tan(0)}^{\tan(\frac{\pi}{2})}\big{(}\dfrac{1}{(\frac{3}{M-1})\bar{\gamma}+(2+(\frac{3}{M-1})\bar{\gamma})t^{2}}\big{)}dt=&\\\
&\int_{\tan(0)}^{\tan(\frac{\pi}{2})}\dfrac{1}{a+bt^{2}}dt\\\ \end{split}$
(83)

To solve the indefinite integral $\int\dfrac{1}{a+bt^{2}}dt$, let
$\tan(\theta)=\frac{\sqrt{b}}{\sqrt{a}}x$, $x=\sqrt{\frac{a}{b}}\tan(\theta)$,
$dx=\sqrt{\frac{a}{b}}\frac{1}{\cos^{2}(\theta)}d\theta$,
$\theta=\tan^{-1}(\sqrt{\frac{b}{a}}x)$,
$\cos(\theta)=\frac{\sqrt{a}}{\sqrt{a+bx^{2}}}$, and
$a+bx^{2}=\frac{a}{\cos(\theta)}$.

$\begin{split}\int\dfrac{1}{a+bt^{2}}dt=&\sqrt{\frac{a}{b}}\int\dfrac{\dfrac{1}{\cos^{2}(\theta)}}{\dfrac{a}{\cos^{2}(\theta)}}d\theta\\\
=&\dfrac{1}{a}\sqrt{\frac{a}{b}}\int d\theta\\\
=&\dfrac{1}{a}\sqrt{\frac{a}{b}}\tan^{-1}(\sqrt{\frac{b}{a}}x).\\\
\end{split}$ (84)

An important assumption has been made as $\sqrt{\frac{b}{a}}\approx 1$, then
$\tan(\tan^{-1}(\theta))=\theta$. By applying the last result to (83), we get

$\begin{split}\int_{\tan(0)}^{\tan(\frac{\pi}{2})}\dfrac{1}{a+bt^{2}}dt=&\dfrac{1}{a}\sqrt{\frac{a}{b}}\dfrac{\pi}{2}\\\
=&\dfrac{1}{(\frac{3}{M-1})\bar{\gamma}}\sqrt{\frac{(\frac{3}{M-1})\bar{\gamma}}{2+(\frac{3}{M-1})\bar{\gamma}}}\dfrac{\pi}{2}\end{split}$
(85) $\begin{split}P_{e-MQAM}(1)=&2(\dfrac{\sqrt{M}-1}{\sqrt{M}})\big{[}1-\\\
&\sqrt{\frac{(\frac{3}{M-1})\bar{\gamma}}{2+(\frac{3}{M-1})\bar{\gamma}}}\big{]}\end{split}$
(86)

Let’s solve the second part of (80) which can reduced in similar to the first
part as follow:

$\begin{split}P_{e-MQAM}(2)=&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\int_{0}^{\frac{\pi}{4}}\big{(}\dfrac{1}{1+(\frac{3}{M-1})\frac{\bar{\gamma}}{2\sin^{2}(\theta)}}\big{)}d\theta\\\
=&(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\int_{0}^{\frac{\pi}{4}}\big{(}\dfrac{(\frac{3}{M-1})\bar{\gamma}}{2\sin^{2}(\theta)+(\frac{3}{M-1})\bar{\gamma}}\big{)}d\theta.\\\
\end{split}$ (87)

let $\tan(\theta)=t$,
$\dfrac{dt}{d\theta}=\sec^{2}(\theta)=1+\tan^{2}(\theta)=1+t^{2}$,
$d\theta=\dfrac{1}{1+t^{2}}dt$, and $\sin(\theta)=\dfrac{t}{\sqrt{1+t^{2}}}$.

$\begin{split}P_{e-MQAM}(2)=&(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}-\\\
&\dfrac{4}{\pi}(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}(\frac{3}{M-1})\bar{\gamma}\\\
&*\int_{\tan(0)}^{\tan(\frac{\pi}{4})}\big{(}\dfrac{1}{(\frac{3}{M-1})\bar{\gamma}+(2+(\frac{3}{M-1})\bar{\gamma})t^{2}}\big{)}dt.\\\
\end{split}$ (88)

Using the solution of the indefinite integral $\int\dfrac{1}{a+bt^{2}}dt$ in
(84) and replace $t=\tan(\theta)$, we get:

$\begin{split}\int_{\tan(0)}^{\tan(\frac{\pi}{4})}\dfrac{1}{a+bt^{2}}dt=&\dfrac{1}{a}\sqrt{\frac{a}{b}}\tan^{-1}(\sqrt{\frac{b}{a}}\tan(\theta))\big{|}_{0}^{\frac{\pi}{4}}\\\
=&\dfrac{1}{(\frac{3}{M-1})\bar{\gamma}}\sqrt{\frac{(\frac{3}{M-1})\bar{\gamma}}{2+(\frac{3}{M-1})\bar{\gamma}}}\\\
&*\tan^{-1}\big{(}\sqrt{\frac{2+(\frac{3}{M-1})\bar{\gamma}}{(\frac{3}{M-1})\bar{\gamma}}}\big{)}\end{split}$
(89)

Therefore,

$\begin{split}P_{e-MQAM}(2)=&(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\big{[}1-\sqrt{\frac{(\frac{3}{M-1})\bar{\gamma}}{2+(\frac{3}{M-1})\bar{\gamma}}}\\\
&*\dfrac{4}{\pi}\tan^{-1}\big{(}\sqrt{\frac{2+(\frac{3}{M-1})\bar{\gamma}}{(\frac{3}{M-1})\bar{\gamma}}}\big{)}\big{]}\\\
=&(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\big{[}1-\sqrt{\frac{1.5\bar{\gamma}}{M-1+1.5\bar{\gamma}}}\\\
&*\dfrac{4}{\pi}\tan^{-1}\big{(}\sqrt{\frac{M-1+1.5\bar{\gamma}}{1.5\bar{\gamma}}}\big{)}\big{]}\\\
\end{split}$ (90)

The total symbol error probability for M-QAM is:

$\begin{split}P_{e-MQAM}=&2(\dfrac{\sqrt{M}-1}{\sqrt{M}})\big{[}1-\sqrt{\frac{1.5\bar{\gamma}}{M-1+1.5\bar{\gamma}}}\big{]}\\\
&(\dfrac{\sqrt{M}-1}{\sqrt{M}})^{2}\big{[}1-\sqrt{\frac{1.5\bar{\gamma}}{M-1+1.5\bar{\gamma}}}\\\
&*\dfrac{4}{\pi}\tan^{-1}\big{(}\sqrt{\frac{M-1+1.5\bar{\gamma}}{1.5\bar{\gamma}}}\big{)}\big{]}\\\
\end{split}$ (91)

##  VIII Conclusion

In this paper, we have presented a unified approach for deriving the
probability of error formulations applicable to BPSK, 16-QAM, and 64-QAM in
Rayleigh fading channels. By considering the statistical properties of
Rayleigh fading and the modulation characteristics of each scheme, our
approach provides a comprehensive framework for analyzing error performance
across different modulation schemes. Importantly, our approach applies to
solving all types of integration processes, producing accurate solutions that
account for the complexities of Rayleigh fading and modulation schemes. By
establishing a unified framework, we have simplified the analysis process and
facilitated a deeper understanding of error behavior, enabling researchers,
engineers, and practitioners to design and optimize communication systems with
confidence. Looking ahead, further research can explore additional modulation
schemes, channel models, and advanced signal processing techniques to extend
the applicability and effectiveness of our approach. The findings presented
herein serve as a foundation for future research and development efforts aimed
at achieving reliable communication in challenging fading environ

##  Appendix A Common Used Functions

The complimentary distribution function (CDF) is calculated by:

$F_{x}(x)=\dfrac{1}{\sqrt{2\pi\sigma^{2}}}\int_{-\infty}^{x}e^{\dfrac{-(\alpha-\mu)^{2}}{2\sigma^{2}}}d\alpha$
(92)

The complimentary error function is:

$\begin{split}erfc\big{(}x\big{)}=&\dfrac{2}{\sqrt{\pi}}\int_{x}^{\infty}e^{-t^{2}}dt\\\
=&1-erf(x),\\\ \end{split}$ (93)

and error $erf(x)$ function is expressed by:

$erf\big{(}x\big{)}=\dfrac{2}{\sqrt{\pi}}\int_{0}^{x}e^{-t^{2}}dt,$ (94)

the derivative of $erf(x)$ function is:

$\dfrac{d}{dx}erf\big{(}x\big{)}=-\dfrac{2}{\sqrt{\pi}}e^{-x^{2}}$ (95)

And also the definite integral of $erfc(x)$, simply calculating it by parts
($\int udv=uv-\int vdu$). By taking $u=erfc(x)$,
$du=\dfrac{2}{\sqrt{\pi}}e^{-t^{2}}$, $dv=dt$, and $v=x$.

$\int erfc\big{(}x\big{)}dx=xerfc(x)-\dfrac{2}{\sqrt{\pi}}\int xe^{-x^{2}}dx,$
(96)

By taking $u=-x^{2}$ and $du=-2xdx$.

$\begin{split}-\dfrac{2}{\sqrt{\pi}}\int
xe^{-x^{2}}dx=&\dfrac{1}{\sqrt{\pi}}\int e^{u}du,\\\
=&\dfrac{1}{\sqrt{\pi}}e^{u}+C\\\ =&\dfrac{1}{\sqrt{\pi}}e^{-x^{2}}+C\\\
\end{split}$ (97) $\int
erfc\big{(}x\big{)}dx=xerfc(x)+\dfrac{1}{\sqrt{\pi}}e^{-x^{2}}+C$ (98)

##  Appendix B The Relation Between The $Q(x)$ and $erfc(x)$ Functions

The Q-function CDF is:

$\begin{split}Q(x)=&Pr[X>x]=1-F(x)\\\
=&\dfrac{1}{\sqrt{2\pi}}\int_{x}^{\infty}e^{\dfrac{-u^{2}}{2}}du\\\
=&\dfrac{1}{2}erfc\big{(}\sqrt{\frac{x}{2}}\big{)}\\\
=&\dfrac{1}{2}[1-erf\big{(}\sqrt{\frac{x}{2}}\big{)}],\\\ \end{split}$ (99)

By taking $z^{2}=\dfrac{u^{2}}{2}$, $z=\dfrac{u}{\sqrt{2}}$,
$dz=\dfrac{1}{\sqrt{2}}du$, and $du=\sqrt{2}dz$. The integration limits will
be as follows: $u=x\Rightarrow z=\dfrac{x}{\sqrt{2}}$ and $u=\infty\Rightarrow
z=\infty$.

$\begin{split}Q(x)=&\dfrac{\sqrt{\pi}}{\sqrt{2\pi}}\int_{\frac{x}{\sqrt{2}}}^{\infty}e^{-z^{2}}dz\\\
=&\dfrac{1}{2}erfc\big{(}\frac{x}{\sqrt{2}}\big{)}\end{split}$ (100)

To get the Q-function from erfc function, we can follow the simple below
calculations:

$\begin{split}erfc\big{(}x\big{)}=&\dfrac{2}{\sqrt{\pi}}\int_{x}^{\infty}e^{-t^{2}}dt\\\
=&1-erf(x),\\\ \end{split}$ (101)

By taking $\dfrac{z^{2}}{2}=u^{2}$, $z=u\sqrt{2}$, $dz=\sqrt{2}du$, and
$du=\dfrac{1}{\sqrt{2}}dz$. The integration limits will be as follows:
$u=x\Rightarrow z=x\sqrt{2}$ and $u=\infty\Rightarrow z=\infty$.

$\begin{split}erfc\big{(}x\big{)}=&\dfrac{2}{\sqrt{\pi}}\int_{x}^{\infty}e^{-u^{2}}du\\\
=&\dfrac{2}{\sqrt{2\pi}}\int_{\sqrt{2}x}^{\infty}e^{-\frac{z^{2}}{2}}dz,\\\
=&2Q(x\sqrt{2}),\\\ \end{split}$ (102)

## Acknowledgment

The authors wish to thank Krishna Sankar who is Principal Systems Design
Engineer of the HMicro, Bangalore area, India for the informative webpage on
wireless communications and simplified explanation of BER of BPSK, PAM, and
QAM in AWGN and Rayleigh fading channels [8].

## References

  * [1] Proakis, _Digital Communications 5th Edition_. McGraw Hill, 2007.
  * [2] J. R. Barry, E. A. Lee, and D. G. Messerschmitt, _Digital Communication_ , 3rd ed. Springer, 2004.
  * [3] E. Säckinger, _Analysis and Design of Transimpedance Amplifiers for Optical Receivers_. John Wiley & Sons, Ltd, 2017, ch. Appendix A: Communication Signals, pp. 397–411.
  * [4] D. Lu, X. Zhou, Y. Yang, J. Huo, J. Yuan, K. Long, C. Yu, A. Lau, and C. Lu, “Theoretical analysis of pam-n and m-qam ber computation with single-sideband signal,” _Science China Information Sciences_ , vol. 64, no. 8, Jun. 2021.
  * [5] M. K. Simon and M.-S. Alouini, _Digital Communication Over Fading Channels_. John Wiley & Sons, Ltd, 2004, ch. 8, pp. 193–258.
  * [6] O. Oyetola, A. Okubanjo, A. A. Okandeji, P. O. Alao, M. O. Osifeko, and O. G. Olasunkanmi, “Symbol error probability of 16-qam system over awgn and rayleigh fading channels,” _African Journal of Science and Nature_ , 2020\.
  * [7] P. Fortune, L. Hanzo, and R. Steele, “On the computation of 16-qam and 64-qam performance in rayleigh-fading channels,” _IEICE Transactions, Japan_ , vol. E75-B, no. 6, pp. 466–475, June 1992.
  * [8] K. Sankar. (2024) Dsp log signal processing. [Accessed: 24 June 2024]. [Online]. Available: https://dsplog.com/