\section{\label{S:Intro}Introduction}
In early 1990's, carbon nanotechnology revolution introduced a plethora of new advanced materials among which the carbon nanotube was notoriously attractive for making nanoscale field effects devices, including vacuum devices. Many labs studied effects associated with field emission from a single CNT or arrays with counted number of isolated CNTs.\cite{Motorola1, Thales, BrightCNT, HotCNT_PRL, CNT_arrays_apps1, CNT_arrays_apps2} Control over fabrication and emission of single CNT field emission devices was excellent and many field emission devices were demonstrated, e.g. field emission radio\cite{CNTradio} or field emission transistor\cite{CNT_FE_transistor}, amplifier\cite{CNT_FE_amp}, and many others.\cite{CNT_arrays_apps1,CNT_arrays_apps2}

\begin{figure}
\includegraphics[width=5.5cm]{Fig_nonuniform}
\caption{Typical micrograph showing large beam transverse spread and nonuniformity. The blue circle marks the cathode's position behind the imaging screen and its size. }\label{F:nonuniform}
\end{figure}

In order to increase the output power, macrosocopic large area CNT fibers, films, yarns and fabrics started to be used to increase the operating currents from pico-/nano- to many amperes. Here, CNTs were thought to replace legacy velvets.\cite{CNTvsVelvet} The multiple benefits of CNT fibers over legacy technology are low turn-on voltage and high emission current at relatively low operating electric field due to inherent high field enhancement factor, and high electrical and thermal conductivities.\cite{fabrication} It was conventionally assumed that emission would be uniform, i.e. uniformity would translate from previously studied arrays of counted CNTs to the large area CNT fibers. However, recent studies that employed field emission microscopy illustrated that emission is never uniform and moreover that the emission area is a function of the electric field (making it cumbersome for calculating current densities.) Fig.\ref{F:nonuniform}, reproduced from our past work,\cite{mypaper} highlights other important issue of the large transverse spread of the emitted beam where beam lands on the imaging screen millimeters away from the physical location of the cathode source (blue circle) after travelling only a millimeter between the cathode and the anode. This clearly points out a very large emittance and therefore very low brightness, making CNT fiber cathodes impractical for applications like rf or microwave traveling wave tubes (operating in GHz range), microscopy and bright X-ray sources for medicine or active scanning. Another issue arises from that---because all the current emerges from a few active spots, it leads to local heating, microbreakdowns\cite{mypaper} and short-lived cathodes.

After experimenting with many fiber arrangements, we found that it is tiny singular fibrils (comprising braided fibers) that set loose due to thermal and field related stress\cite{mypaper,morphology_dependent} and that eventually focus the field due to their high aspect ratio and become point-like randomized intense electron emitters eventually exploding (seen as micro-breakdowns) and re-populating surrounding areas with more new-born fibrils. This process repeats itself until the cathodes stops operating while emission always look like a family of single electron rays going in many directions that are not aligned with the desired main longitudinal propagation direction, such as in Fig.\ref{F:nonuniform}. To mitigate this issue, hypothesised to be the major problem, we study new cathode production technology where fibers are electroplated with Ni and laser cut; all to suppress the fibril occurrence and regeneration. Through experimental measurements and electrostatic and beam dynamics modeling, emission uniformity and beam brightness were analyzed.

\section{\label{S:experimental}Experimental}
To prepare the field emission CNT fiber cathodes, a commercially available CNTs fiber from DexMat, Inc was used. The fiber is made by a wet-spinning technology \cite{fabrication}---pre-grown arrays of CNTs are dissolved in an acid to form a spinnable liquid dope that is extruded through a spinneret into coagulant bath to remove acid, and then dried in an oven. The resulting product is highly aligned and densely packaged CNTs in a form of a fiber. DexMat fibers have high electrical and thermal conductivity. Such fibers were shown to feature anisotropic field emission\cite{unisotropic_emission} that is emission takes place along the fiber (not from side walls) which is a great property allowing for control over emittance. Raman spectroscopy shows the G peak positions at 1583 cm$^{-1}$ suggesting rich crystalline graphitic content as expected from high quality fibers (Fig.\ref{F:sample}A).

To mitigate the described stray fibril problem, few fibers of the described kind were placed side by side and electroplated with Ni in an electrochemical bath and flush cut from the top to the required length of about $5\,\text{mm}$ with a femtosecond micromachining laser beam. Then, it was welded on $1\times 1\,\text{inch}$ Ni base. The final fabricated structure can be seen in Fig.\ref{F:sample}B.

We tested two samples referred to as Sample A and Sample B through the rest of the paper. In Sample A fibers were twisted and in Sample B were not, i.e. simply place along each others. As scanning electron microscopy (SEM) demonstrates in Fig.\ref{F:sample}C, additional fiber twisting enabled a dense core in Sample A, while Sample B (Fig.\ref{F:sample}D) has visible voids between the individual fibers. Otherwise, both sample have fiber core diameter of $\sim 150\,\mu\text{m}$ and Ni shell thickness of $\sim50\,\mu\text{m}$. Sample A and Sample B have height of $4.8\,\text{mm}$ and $4.6\,\text{mm}$ respectively (see Fig.\ref{F:sample}B).

\begin{figure}
    \centering
    \includegraphics[scale=0.045]{Fig_sample}
    \caption{(A) Raman spectra of the cathode surface showing a crystalline graphitic peak. (B) Electroplated CNT fiber welded on a Ni base; $\text{H}=4.8\,\text{mm}$ for Sample A and $\text{H}=4.6\,\text{mm}$ for Sample B. (C) and (D) SEM images of Sample A and Sample B, where scale bars are 50 $\mu$m.}
    \label{F:sample}
\end{figure}

DC current tests and field emission microscopy were performed in our custom field emission microscope described in great detail in Ref.\onlinecite{exp_setup}. Images were processed by a custom image processing algorithm FEpic described elsewhere.\cite{mypaper2}

\section{\label{S:ivrelation}Field Emission Imaging and Conditioning}

After sample was installed and gap was tuned using a doublet of two orthogonal optical microscopes, the physical location of the fiber is determined and labeled. To do that, the test chamber is illuminated. Because the imaging anode YAG screen is semitransparent, the location of the fiber can be immediately seen and captured by photographing. The core of the fiber is marked with red circle for reference in Figs.\ref{F:sampleA}A and \ref{F:sampleB}A.

After that voltage is applied and field emission images are taken concurrently with I-V curves. Fig.\ref{F:sampleA}B shows emission micrograph of the Sample A. The improvement is immediately obvious when compared with Fig.\ref{F:nonuniform}. First, the emission spot appears exactly at the optical projection of the cathode. This means beam divergence angle is small, so emittance can be expected to be low. Second, there is only a single spot and its size is comparable to the size of the fiber core---this is an indication of uniformity and small angular spread of the electron beam. The same exact behaviour was observed for Sample B, as given in Fig.\ref{F:sampleB}. No evidence suggesting the stray fibril issue was observed for neither cathode. These results highlight that such a simple electroplating strategy is extremely effective at yielding emission uniformity and spatial coherence, thereby boosting the transverse beam brightness.

\begin{figure}
\includegraphics[scale=0.47]{Fig_sampleA}
\caption{A) Sample A seen through the YAG screen when the lights is on in the chamber. Its fiber core is marked with red circle. B) FE micrograph of the same region at the gap of $200\,\mu\text{m}$.}\label{F:sampleA}
\end{figure}

\begin{figure}
\includegraphics[scale=0.47]{Fig_sampleB}
\caption{A) Sample B seen through the YAG screen when the lights is on in the chamber. Its fiber core is marked with red circle. B) FE micrograph of the same region at the gap of $200\,\mu\text{m}$.}\label{F:sampleB}
\end{figure}

\begin{figure}
\includegraphics[width=7.cm]{Fig_Condition}
\caption{A) The conditioning scheme of Sample A. Both ramp up and down curves are shown. There is a clear decrease in performance. B) Conditioning scheme of Sample B. Only ramp-up curves are shown. There is no considerable change in performance.}\label{F:condition}
\end{figure}

Fig.\ref{F:condition} shows cathode conditioning I-V curves for Sample A and B at interelectrode gap of 200 $\mu$m. Conditioning\cite{tungsten_conditioning,diamond_conditioning,mypaper} or cycling, where the applied voltage is ramped up and down to a progressively higher number in every consecutive cycle until the desired operating current is achieved, is a crucial procedure to maximize field emission cathode performance and ensure stability at the operating point.

We found that the electroplated fibers should be conditioned with small incremental steps to avoid adverse effects such as sudden burn-down. Fig.\ref{F:condition}A demonstrates the case where the maximal field was doubled with respect to the previous conditioning cycle, such that the current went up from 3 to 22 $\mu$A. Next, after completing the ramp down Sample A stopped working completely, which was possibly due to applying electric power that exceed that in the previous run by more than an order of magnitude. Its operation could not be rejuvenated by applying higher electric fields. This is unlike a conditioning scheme that was used for Cathode B as shown in Fig.\ref{F:condition}B. The emission current was doubled at every conditioning cycle: up to 1 $\mu$A and down to 0, then to 2 $\mu$A, to 4 $\mu$A, to 8 $\mu$A, to 16 $\mu$A, and finally to 32 $\mu$A. By doing so, Cathode B was conditioned $softly$ (compared to Sample A) maintaining and enhancing its performance: $i$) the resulting operating field went up and doubled, reaching same exact value where Cathode A burned down; $ii$) turn-on field and field enhancement factor remained nearly the same meaning that Cathode B was conditioned to stably sustain higher local field.

When it is compared to our past cathode designs, detailed in Ref.\onlinecite{mypaper}, they emit less at any given field. This is an expected result because (with stray fibrils mitigated) the field enhancement is reduced. However, turn-on fields are still very low, between 1 and 2.5 V/$\mu$m. Because the beam was tight suggesting high current density we limited our measurements to between 10--100 $\mu$A as the power density deposition at the imaging screen could attain above 1 kW/cm$^2$ at the voltage source limit of 1100 V, thereby literally drilling holes in it.\cite{power-density-1, power-density-2} At 1100 V, Sample A maxed out at 20 $\mu$A and Sample B 30 $\mu$A, respectively. Again both cathodes had similar metrics. Having these metrics and qualitative results in mind, a step was taken to carry out more quantitative analysis and calculate cathodes' emittance and brightness. All detailed in the next section.

\section{\label{S:emittance}Emittance and Brightness}

In the phase space $(x,x')$, $x$ is spatial position and $x'=\frac{dx}{dz}=\frac{dx/dt}{dz/dt}=\frac{v_x}{v_z}$ is the slope of the trajectory from longitudinal centrosymmetric axis of each particle. Then, rms emittance $\tilde{\epsilon}_x$ is defined as
\begin{equation}\label{E:rmsemit}
    \tilde{\epsilon}_x = \sqrt{\langle \Delta x^2 \rangle \langle \Delta x'^2 \rangle - \langle \Delta x \Delta x' \rangle}
\end{equation}
where $\Delta x=x-\langle x\rangle$ and $\Delta x'=x'-\langle x'\rangle$. For a beam with cylindrical symmetry in $(x,y)$ and $(x',y')$ centered around zero, $\langle x\rangle$ and $\langle x'\rangle$ are zero. The cross-corelation term, $\langle \Delta x \Delta x' \rangle$ can be removed with proper beam optics.\cite{jarvis_thesis} Then, Eq.\ref{E:rmsemit} reduces to
\begin{equation}\label{E:rmsemit2}
    \tilde{\epsilon}_x = \sqrt{\langle x^2 \rangle \langle x'^2 \rangle}=\sigma_x\sqrt{\langle x'^2 \rangle}.
\end{equation}

Rms emittance is a function of the beam energy as $x'$ is changing under acceleration, and is not useful while comparing beam or beams at different energies. On the other hand, from Liouville's theorem, normalized emittance is a conserved quantity under acceleration as long as the beam is only subjected to conservative forces. The relation between rms emittance and normalized emittance is given by
\begin{equation}\label{E:normemit}
    \epsilon_x^\text{N} = \gamma\beta\tilde{\epsilon}_x
\end{equation}
where $\gamma = \frac{1}{\sqrt{1-\beta^2}}$ is the Lorentz factor and $\beta = \frac{v}{\text{c}} \approx \frac{v_z}{\text{c}}$. In our case, $\gamma\approx1$ because energy is $<=$1 keV. Mean-transverse energy, $\text{MTE}$, is $\frac{1}{2}\text{m}_\text{e} \langle v^2 \rangle$, where $v^2=v_x^2+v_y^2$. Because of the cylindrical symmetry in $(v_x,v_y)$, $\text{MTE}\approx\frac{1}{2}\text{m}_\text{e} \langle 2v_x^2 \rangle=\text{m}_\text{e}\langle v_x^2\rangle$. Then, after substituting Eq.\ref{E:rmsemit2}, in terms of $\text{MTE}$, Eq.\ref{E:normemit} becomes
\begin{equation}\label{E:fundametal}
    \epsilon_x^\text{N} = \frac{v_z}{\text{c}}\cdot\sigma_x\cdot\sqrt{\langle \frac{v_x^2}{v_z^2} \rangle}=\sigma_x\cdot\sqrt{\frac{\text{MTE}}{\text{m}_\text{e}\text{c}^2}}
\end{equation}
Practically, the normalized emittance at the cathode surface is calculated as follows. If the radius of the uniformly emitting surface of the cathode is $r_\text{i}$, $\sigma_x\approx r_\text{i}$ can be taken. Moreover, at the surface, MTE is due to statistical distribution of electrons inside the cathode itself. So, it is intrinsic, and is further redefined as $\text{MTE}_\text{i}$. Then, Eq.\ref{E:fundametal} becomes
\begin{equation}\label{E:ex_finalform}
    \epsilon_x^\text{N}=r_\text{i}\cdot\sqrt{\frac{\text{MTE}_\text{i}}{\text{m}_\text{e}\text{c}^2}}
\end{equation}

By using normalized emittance, normalized transverse brightness, $B_\text{N}$, can be calculated as
\begin{equation}\label{E:brightness}
 B_\text{N}= \frac{2\,I}{\epsilon_x^\text{N} \epsilon_y^\text{N}}
\end{equation}
where $I$ is the emitted current. $\epsilon_x^\text{N}=\epsilon_y^\text{N}$ can be taken in cylindrical symmetry.

\begin{figure}
\includegraphics[scale=0.24]{Fig_moveaway}
\caption{Series of field emission micrographs of Sample B with the screen moved away progressively from $200\,\mu\text{m}$ gap to $1600\,\mu\text{m}$ gap with $200\,\mu\text{m}$ steps. At each step, a micrograph is captured.}\label{F:moveaway}
\end{figure}

To determine $r_\text{i}$ and $\text{MTE}_\text{i}$ in Eq.\ref{E:ex_finalform}, we conducted experimental measurements in combination with beam dynamics in GPT (General Particle Tracer).\cite{gpt_ref} In the measurements, the imaging screen is moved away from the cathode progressively. The voltage is set accordingly to keep the current constant at 20 nA to enable strong beam image signal but avoid additional beam expansion due to vacuum space charge effect. A micrograph at each step is recorded  (see Fig.\ref{F:moveaway}). Increase in the spot size due to larger time of flight is measured. As it is seen in Fig.\ref{F:moveaway}, the spots are Gaussian in nature with dense center and faint tails. Each spot can be modeled mathematically with a cylindrically symmetric Gaussian as
\begin{equation}\label{E:gauss_fit}
    p=A\cdot \text{exp}\left( -\frac{(x-x_c)^2+(y-y_c)^2}{2\sigma_{\text{spot}}^2}\right)+C
\end{equation}
to extract projected transverse beam size. Here, $A$ is the amplitude, $\sigma_\text{spot}$ is the standard deviation, $C$ is the background offset, $p$ is the intensity, $(x,y)$ are the space dimensions, $(x_c,y_c)$ are the coordinates of the peak.\cite{mypaper2} The model parameters $A$, $\sigma_\text{spot}$, and $C$ for each spot are computed with least-square fitting method. After fitting, the emission spot diameter is taken as $2\sigma_\text{spot}$. An exemplary 3D fitting done by FEpic for the beam imaged at 200 $\mu$m gap is presented in Fig.\ref{F:gauss}; with the black mesh surface being the fitting surface. The resulting dependence of the spot size, first measured (Fig.\ref{F:moveaway}) and processed by FEpic, versus distance is shown in Fig.\ref{F:gapvsize} with black solid circles. The data in the figure is only for Sample B. Because Sample A burned down, studies similar to those presented in Fig.\ref{F:moveaway} could not be carried out.

\begin{figure}
\includegraphics[scale=0.2]{Fig_gauss}
\caption{The color surface shows the beam spot for $200\,\mu\text{m}$ gap in 3D. The black mesh surface shows its mathematical fit in Eq.\ref{E:gauss_fit}.}\label{F:gauss}
\end{figure}

\begin{figure}
\includegraphics[width=6.5cm]{Fig_gapvsize}
\caption{Comparison of experimental and computational final beam spot size as the screen is moving away from the cathode. $\text{MTE}_\text{i}$ of $250\,\text{eV}$ and $r_\text{i}$ of $75\,\mu\text{m}$ were used in GPT modeling.} \label{F:gapvsize}
\end{figure}

\begin{figure}
\includegraphics[width=6.5cm]{Fig_comsol}
\caption{Illustration of electric field computed in COMSOL for a $1\,\text{mm}$ gap. The color plot shows the field magnitude and contour. The dark region is the fiber core, and the gray region is the Ni shell.}\label{F:comsol}
\end{figure}

\begin{figure}
\includegraphics[width=8.5cm]{Fig_gpt}
\caption{In GPT: (A) Initial uniform beam distribution at the cathode surface in real-space, where $r_\text{i}$ is the radius of the beam. (B) The initial distribution in momentum-space, where $\beta_x=v_x/c$, $\beta_y=v_y/c$, and $r_\beta$ is the radius. Final distribution in (C) real space and (D) momentum space when the screen at is $1\,\text{mm}$.}\label{F:gpt}
\end{figure}

To calculate the phase space volume for emittance and brightness estimations and obtain $r_\text{i}$ and $\text{MTE}_\text{i}$, we switched to beam dynamics in GPT by comparing it with the experiment in Fig.\ref{F:moveaway}. To do that, a field map for each interelectrode gap was computed in COMSOL by solving Poisson's equation with given boundary conditions. An exemplary field distribution for 1 mm gap is shown in Fig.\ref{F:comsol}. Then, the field maps were imported to GPT. In GPT, the initial particle distribution in the real space $(x,y)$ (Fig.\ref{F:gpt}A)  and momentum space $(\beta_x,\beta_y)=(v_x/\text{c},v_y/\text{c})$ (Fig.\ref{F:gpt}B) has to be set to propagate the beam. Because the cathode itself circular, we used a circular uniform distribution with the radius $r_\text{i}$ in position-space and the radius $r_\beta$ in $\beta$-space. $\text{MTE}_\text{i}$ of the fiber is expected to be $250\,\text{meV}$.\cite{mte_fairchild} To be used in GPT, $\text{MTE}_\text{i}$ is converted into $r_\beta$ as
\begin{equation}\label{E:mte_vs_maxB}
    r_\beta=\sqrt{\frac{4\,\text{MTE}_\text{i}}{\text{m}_\text{e}\text{c}^2}}.
\end{equation}
$\text{MTE}_\text{i}$ of 250 meV translates to $r_\beta$ of $1.4\times10^{-3}$. At the same time, $r_\text{i}$ remains a free model parameter to be found by finding the best agreement between GPT with the experiment.

The beam was launched at a charge corresponding to 20 nA and allowed to drift through the distance corresponding to a specified interelectrode gap. Here, the beam dynamics is computed self-consistently taking the COMSOL calculated field. The final distribution in $(x,y)$ and $(\beta_x,\beta_y)$ was captured at a distance corresponding to the imaging YAG screen of the microscope and are shown in Fig.\ref{F:gpt}C and D. In the $(x,y)$ space, the standard deviation $\sigma_{gpt}$ of such projections were calculated for every cathode-anode gap, and $2\sigma_{gpt}$ was taken as the computed beam diameter (analogous to FEpic image processing of the experimental images). It was established that in GPT, when $r_\text{i}$ was set to the physical radius of the fiber core of $75\,\mu\text{m}$, and $r_\beta$ was set by MTE of 250 meV,\cite{mte_fairchild} then the final diameter (spot size $2\sigma_{gpt}$) of the resulting beam projection was in a very good quantitative agreement with the experiments (see Fig.\ref{F:gapvsize}). Note, the GPT results were fairly insensitive to MTE values set between 25 (typical Fermi level value for CNT) and 250 meV, and magnification was due to the radial field distribution. This points out that the divergence between the experiment and GPT (seen for the gap ranging between 1 and 1.6 mm) stems from the difference between the idealized computed and the actual field distribution in the gap. The summary of the results in Fig.\ref{F:gapvsize} confirms that the entire whole fiber surface actively and uniformly emitting with a small $\text{MTE}$.

Finally, using Eq.\ref{E:ex_finalform} and substituting $r_\text{i}=75\,\mu\text{m}$ and $\text{MTE}_\text{i}=250\,\text{eV}$, the upper limit of the normalized emittance on the fiber cathode surface can be estimated as
\begin{equation}\label{E:num_emit}
    \epsilon_x^N=0.052\,\text{mm}\,\text{mrad} = 52\,\text{nm}.
\end{equation}
From this, taking the measured current (limited to 10--100 $\mu$A due to extremely high power density), as shown in Fig.\ref{F:condition}, the normalized brightness for 50 $\mu$A dc current is $B_N=3.7\times10^{10}\frac{\text{A}}{\text{m}^2\text{rad}^2}$. The same very fibers can draw currents of 1--10 A when operated in pulsed mode with a pulse length of 100--300 ns.\cite{fiber-IVEC2020} Using the estimated emittance of $52\,\text{nm\,rad}$, the brightness in the pulsed mode, preferable mode in most VED HPM applications, attains a notable value of $B_N=4.4\times10^{15}\frac{\text{A}}{\text{m}^2\text{rad}^2}$. This number is outstanding and is comparable with brightness metrics in the the state-of-the-art microwave/rf accelerator injectors.\cite{SOTA_Br} This (phase space) brightness can be converted into geometrical brightness, a definition of brightness commonly employed in the electron microscopy literature. The geometrical reduced brightness is defined as $B_r^G=\frac{\text{d}I}{\text{d}\Omega}\frac{1}{U}\frac{1}{S_{\text{cathode}}}$, where $\Omega$ is the solid angle, $U$ is the voltage at which the current $I$ is measured, and $S_{\text{cathode}}$ is the emission area of the cathode. Our calculations show that in pulsed mode it could attain $B_r^G=5.7\times 10^7$ A m$^{-2}$sr$^{-1}$V$^{-1}$. This number is within the range obtained for single CNT emitters.\cite{BrG_PRl2005}

\section{\label{S:conclusion}Conclusion}

In conclusion, we presented a simple and efficient field emission cathode design where CNT fiber core was plated with nickel shell. This design had two important functions. First, it compresses the core, provides mechanical strength thereby preventing stray fibril formation during conditioning and operation. Second, such design (while slightly reducing field enhancement and increasing turn-on field) reduces the fringing field on the CNT fiber and therefore the defocusing radial field.

As field emission microscopy directly demonstrated, both tested cathodes featured excellent spatially coherent emission. Field emission microscopy aided by image processing and beam dynamics simulations confirmed that the entire fiber core of 150 $\mu$m in diameter actively and uniformly emitted electrons, as well as enabled phase space analysis. All of these combined allowed to quantify the observed emission coherence through calculating emittance and brightness. The extremely low emittance resulting in record brightness highlight a simple and practical path forward for the CNT fiber technology that has long been expected to advance high frequency vacuum power devices but had limited success due to low brightness.

Finally and most importantly, it was demonstrated that the nanoscopic single CNT cathode technology can be translated to the macroscopic fiber CNT level in terms of emission uniformity. In other words, spatial coherence and uniformity (intrinsic to a single CNT emitter) can be achieved in a CNT fiber comprised out of billions of single CNT's. The obtained brightness figures of merit further confirm this technology translation in that ultimate single CNT emitter brightness is feasible to attain for CNT fiber cathodes.


\section{Introduction}
It is well known that the classical coding theory studies perfect codes in Hamming graphs and these graphs are distance-transitive. In $1973$, Biggs \cite{Biggs} initiated an investigation of perfect codes in distance-transitive graphs. Since distance-transitive graphs are a family of distance-regular graphs, it is reasonable to study perfect codes in distance-regular graphs. Neumaier \cite{N} introduced the notion of a completely regular code and proved that a perfect code in a distance-regular graph is indeed a completely regular code. We refer to the monograph \cite[Chap.~11]{BCN} for more background on perfect codes in distance-regular graphs. In this paper, we study the problem that which of distance-regular graphs admit a perfect $1$-code. In some literature, an efficient domination set is used instead of a perfect $1$-code (see for example \cite{CMO} and \cite{CMOS}). For an overview on recent progress of this topic, we refer to \cite{CMOS}. We first state some observations and equations on perfect $1$-codes and then we characterize distance-regular line graphs which admit a perfect $1$-code. Furthermore, we state some facts about perfect codes in antipodal distance-regular graphs and give an overview on perfect $1$-codes in distance-regular graphs with small diameter at most $4$. Moreover, we characterize all known distance-regular graphs with small valency at most $4$, the distance-regular graphs with known putative intersection arrays for valency $5$, and all distance-regular graphs with girth $3$ and valency $6$ or $7$ which admit a perfect $1$-code.
\section{Preliminaries}
In this paper, all graphs are undirected and simple, i.e., there are no loops or multiple edges. Moreover, we consider the eigenvalues of the adjacency matrix of a graph. A connected graph $\Gamma$ is called distance-regular with diameter $d$ and intersection array
\begin{equation*}
\{b_{0},b_{1},\ldots,b_{d-1};c_{1},c_{2},\ldots,c_{d}\}
\end{equation*}
whenever for each pair of vertices $x$ and $y$ at distance $i$, where $0 \leq i \leq d$, the number of neighbours of $x$ at distance $i+1$ and $i-1$ from $y$ are constant numbers $b_{i}$ and $c_{i}$, respectively. This implies that a distance-regular graph is regular with valency $b_{0}=k$ and the number of neighbours of $x$ at distance $i$ from $y$ is a constant number $k-b_{i}-c_{i}$ which is denoted by $a_{i}$. A $k$-regular graph with $n$ vertices is called strongly regular with parameters $(n,k,\lambda,\mu)$ whenever the number of common neighbours of two adjacent vertices is $\lambda$ and the number of common neighbours of two non-adjacent vertices is $\mu$. Note that for a strongly regular graph to be of diameter $2$ and thus a distance-regular graph, it needs to be connected and non-complete. Moreover, for a distance-regular graph with diameter $d$, the number of vertices at distance $i$ from an arbitrary given vertex is constant and denoted by $K_{i}$. Furthermore, $$K_{i+1}=\frac{K_{i}b_{i}}{c_{i+1}},$$ where $i=0,1,\ldots,d-1$ and $K_{0}=1$.

Recall that a projective plane of order $q$ is a point-line incidence structure such that each line has $q+1$ points, each point is on $q+1$ lines, and every pair of points is on a unique line. Furthermore, the incidence graph of a projective plane is a bipartite distance-regular graph with diameter three and intersection array $\{q+1,q,q;1,1,q+1\}$. Moreover, the distinct eigenvalues (of the adjacency matrix) of this graph are $\{\pm(q+1),\pm \sqrt{q}\}$.

Let $\Gamma$ be a graph with vertex set $V$. Then any subset $C$ of $V$ is called a code in $\Gamma$. Let $\overline{\Gamma_{t}(c)}$ denote the set of vertices at distance at most $t$ from $c$, where $c \in C$. Then the code $C$ in $\Gamma$ is called perfect $t$-code whenever $\{\overline{\Gamma_{t}(c)} \mid c \in C\}$ is a partition of the vertex set $V$ (cf. \cite{K}). This implies that a code $C$ is perfect $1$-code whenever $C$ is an independent set and every vertex outside $C$ has a unique neighbour in $C$. Furthermore, it is trivial to see that if $C_{1}$ and $C_{2}$ are two perfect $1$-codes in a graph, then $|C_{1}|=|C_{2}|$ since there exists a bijection between $C_{1}$ and $C_{2}$ by the definition of a perfect $1$-code. Moreover, if $C$ is a perfect $1$-code in a $k$-regular graph $\Gamma$, then
\begin{equation} \label{equation}
|C|=\frac{|V|}{k+1},
\end{equation}
because $\{\Gamma_{1}(c)\mid c \in C\}$ is a partition for the vertex set $V$ and each part has size $k+1$.
\begin{rem} \label{completegraph}
Let $\Gamma$ be a regular graph with vertex set $V$. Then $\Gamma$ admits an one-element subset of $V(\Gamma)$ as a perfect $1$-code if and only if $\Gamma$ is a complete graph.
\end{rem}
The following two observations are trivial by the definition of a perfect $1$-code but are useful.
\begin{Observation} \label{Observation1}
Let $C$ be a perfect $1$-code with at least two elements in a connected graph and $x,y \in C$. Then $d(x,y) \geq 3$. Moreover, there exist at least two elements at distance $3$ in $C$. To see this let $x \in C$ and $y$ be a vertex at distance $2$ from $x$ which is indeed outside $C$. Then there exists a unique element $z \in C$ which is adjacent to $y$ and therefore the distance between $x$ and $z$ is $3$.
\end{Observation}
\begin{Observation} \label{Observation2}
Let $C$ be a perfect $1$-code in a connected regular graph with vertex set $V$ and valency $k$. Then $\{C,V \backslash C\}$ is an equitable partition with the quotient matrix $$\begin{bmatrix}0&k\\1&k-1\end{bmatrix}.$$ Therefore $-1$ must be an eigenvalue of (the adjacency matrix of) this graph.
\end{Observation}
\subsection{Completely regular codes}
Let $\Gamma$ be a connected regular graph with vertex set $V$ and a code $C$, where $|C|\geq 2$. Then the number
\begin{equation*}
d(C):=\min\{d(x,y)\mid x,y \in C,x\neq y\}
\end{equation*}
is called the minimum distance of $C$. The distance $v \in V$ from $C$ is defined by
\begin{equation*}
d(v,C):=\min\{d(v,w)\mid w \in C\}
\end{equation*}
and the number
\begin{equation*}
t(C):=\max\{d(v,C)\mid v \in V\}
\end{equation*}
is called the covering radius of $C$. Let
\begin{equation*}
C_{\ell}:=\{v \in V\mid d(v,C)=\ell\},
\end{equation*}
where $\ell=0,1,\ldots,t(C)$. Then the code $C$ is called completely regular whenever for all $\ell$, every vertex in $C_{\ell}$ has the same number $c_{\ell}$ of neighbours in $C_{\ell-1}$, the same number $b_{\ell}$ of neighbours in $C_{\ell+1}$ and the same number $a_{\ell}$ of neighbours in $C_{\ell}$.  It is trivial to see that every one-element code is completely regular in a distance-regular graph. This definition was fist introduced by Neumaier \cite{N}. He proved that a code $C$ in a distance-regular graph is a perfect code if and only if it is a completely regular code with $d(C)=2t(C)+1$ (cf. \cite[Thm.~4.3]{N}). It follows that if $C$ is a perfect $1$-code in a distance-regular graph, then the code $C$ is a completely regular code with $d(C)=3$ since $t(C)=1$.
\section{Antipodal distance-regular graphs} \label{secantipodal}
Let $\Gamma$ be an antipodal distance-regular graph with diameter $d\geq 3$. Then the folded graph of $\Gamma$ which is denoted by $\overline{\Gamma}$ is a graph whose vertex set is the fibers and two fibers are adjacent whenever there exists an edge between them in the graph $\Gamma$. Recall that two fibers have the same size and if two are adjacent in $\overline{\Gamma}$, then there exists a perfect matching between them in the graph $\Gamma$. Moreover, each vertex in one fiber is adjacent to at most one vertex in another fiber since each pair of vertices in a fiber is at distance $d\geq 3$. Let $C$ be a perfect code in $\Gamma$. Then $C$ is a disjoint union of some fibers (cf. \cite[last Remark on p.~349]{BCN}). It follows that if $\overline{C}$ is the corresponding code in the folded graph $\overline{\Gamma}$, then $\overline{C}$ is a perfect code in $\overline{\Gamma}$. Therefore we have the following proposition.
\begin{prop} \label{antipodal}
Let $\Gamma$ be an antipodal distance-regular graph with diameter $d\geq 3$. Then the graph $\Gamma$ admits a perfect code $C$ if and only if $C$ is a disjoint union of some fibers such that $\overline{C}$ is a perfect code in the folded graph $\overline{\Gamma}$.
\end{prop}
Recall that the folded graph of an antipodal distance-regular graph with diameter $3$ is a complete graph. Therefore we can conclude the following corollary about a perfect $1$-code in such a graph.
\begin{cor} \label{antipodaldrg3}
Let $\Gamma$ be an antipodal distance-regular graph with diameter $3$. Then a code $C$ is $1$-perfect if and only if $C$ is a fiber.
\end{cor}
Furthermore, the folded graph of an antipodal distance-regular graph with diameter $d=4$ or $5$ is a strongly regular graph and therefore we can conclude the following corollary about a perfect $1$-code in these graphs since there is no perfect $1$-code in a strongly regular graph.
\begin{cor} \label{antipodaldrg45}
There is no perfect $1$-code in an antipodal distance-regular graph with diameter $d=4$ or $5$.
\end{cor}
We note that the folded graph of the Doubled odd graph $\operatorname{DO_{n}}$ with diameter $2n-1$ is the Odd graph $\operatorname{O_{n}}$ with diameter $n-1$ and therefore we can conclude the following corollary.
\begin{cor} \label{DoubleOdd}
The Doubled odd graph $\operatorname{DO_{n}}$ admits a perfect $1$-code if and only if the Odd graph $\operatorname{O_{n}}$ admits a perfect $1$-code.
\end{cor}
\section{Distance-regular line graphs}
In this section, we characterize distance-regular line graphs which admit a perfect $1$-code. We denote the line graph of a graph $\Gamma$ by $\operatorname{L(\Gamma)}$. The main theorem of this section is as follows.
\begin{thm}
Let $\Gamma$ be a distance-regular graph with least eigenvalue $-2$. Then $\Gamma$ admits a perfect $1$-code if and only if $\Gamma$ is isomorphic to one of the following graphs.
\begin{itemize}
  \item The cycle graph $C_{6n}$,
  \item the line graph of the Petersen graph,
  \item the line graph of the Tutte-Coxeter graph.
\end{itemize}
\end{thm}
We first observe that if $\Gamma$ is a $k$-regular graph with $n$ vertices and the line graph $\operatorname{L(\Gamma)}$ admits a perfect $1$-code $C$, then by Equation \ref{equation}, we have
\begin{equation*}
|C|=\frac{\frac{nk}{2}}{2k-1}.
\end{equation*}
Furthermore, every vertex in $C$ is related to an edge in the graph $\Gamma$ and therefore the perfect $1$-code $C$ can be considered as a $1$-regular induced subgraph, say $\overline{C}$, which is also a vertex cover of the graph $\Gamma$. In other words, the line graph $\operatorname{L(\Gamma)}$ contains $|\overline{C}|$ edges with mutually disjoint closed edge neighborhoods. Recall that the closed edge neighborhood of an edge $e$ consists
of the neighborhood of $e$ together with the edge $e$ itself.  Moreover, $-(k-1)$ must be an eigenvalue of the graph $\Gamma$ because $-1$ must be an eigenvalue of the line graph $\operatorname{L(\Gamma)}$ (cf. Obs. \ref{Observation2}).
Therefore we can conclude the following proposition.
\begin{prop} \label{bipartite}
Let $\Gamma$ be a $k$-regular graph such that its line graph $\operatorname{L(\Gamma)}$ admits a perfect $1$-code. Then $-(k-1)$ is an eigenvalue of the graph $\Gamma$. Moreover, if the graph $\Gamma$ is bipartite, then $\pm k$ and $\pm(k-1)$ must be the eigenvalues of the graph $\Gamma$.
\end{prop}
Distance-regular graphs with least eigenvalue $-2$ have been classified as follows.
\begin{thm}\label{thm:drg-2} \cite[Thm.~3.12.4 and 4.2.16]{BCN} Let $\Gamma$ be a distance-regular graph with least eigenvalue $-2$. Then $\Gamma$ is a cycle of even length, or its diameter $d$ equals $2,3,4,$ or $6$. Moreover,
\begin{itemize}
\item If $d=2$, then $\Gamma$ is a cocktail party graph, a triangular graph, a lattice graph, the Petersen graph, the Clebsch graph, the Shrikhande graph, the Schl\"{a}fli graph, or one of the three Chang graphs,
   \item If $d=3$, then $\Gamma$ is the line graph of the Petersen graph, the line graph of the Hoffman-Singleton graph, the line graph of a strongly regular graph with parameters $(3250,57,0,1)$, or the line graph of the incidence graph of a projective plane,
   \item If $d=4$, then $\Gamma$ is the line graph of the incidence graph of a generalized quadrangle of order $(q,q)$,
   \item If $d=6$, then $\Gamma$ is the line graph of the incidence graph of a generalized hexagon of order $(q,q)$.
   \end{itemize}
\end{thm}
We note that the distance-regular graphs with least eigenvalue {\em larger} than $-2$ are also known. Besides the complete graphs (with least eigenvalue $-1$), there are the cycles of odd length. Recall that for a complete graph, perfect $1$-codes are only one-element subsets of the vertex set (see also Rem. \ref{completegraph}). Furthermore, a cycle graph $C_{n}$ of length $n$ has eigenvalues $2\cos(\frac{2\pi j}{n})$, where $j=0,1,\ldots,n-1$. If this graph admits a perfect $1$-code, then $-1$ must be an eigenvalue of this graph by Observation \ref{Observation2}. Therefore we can conclude the following straightforward proposition.
\begin{prop} \label{Cycle}
A cycle graph $C_{n}$ admits a perfect $1$-code if and only if $3$ divides $n$.
\end{prop}
Now we have to investigate four cases $d=2$, $d=3$, $d=4$ and $d=6$ in Theorem \ref{thm:drg-2}. The case $d=2$ can be ruled out by Observation \ref{Observation1}. For $d=3$, the line graph of the Petersen graph is an antipodal distance-regular graph and every fiber of this graph is a perfect $1$-code (cf. Cor. \ref{antipodaldrg3}). Furthermore, the line graph of the Hoffman-Singleton graph, the line graph of a strongly regular graph with parameters $(3250,57,0,1)$, and the line graph of the incidence graph of a projective plane don't admit a perfect $1$-code because $-1$ is not an eigenvalue of these graphs.
\subsection{The line graph of the incidence graph of a generalized quadrangle}
Let $\Gamma$ be the incidence graph of a generalized quadrangle of order $(q,q)$. Then it has intersection array $\{q+1,q,q,q;1,1,1,q+1\}$ and five distinct eigenvalues $\{\pm(q+1),\pm \sqrt{2q},0\}$ (cf. \cite[Sec.~$6.5$]{BCN}). Therefore, by Proposition \ref{bipartite}, the only possible case is $q=2$. If $q=2$, then the graph $\Gamma$ is indeed the Tutte-Coxeter graph. If the line graph of this graph admits a perfect $1$-code $C$, then $|C|=9$ by Equation \ref{equation}, because it is $4$-regular with $45$ vertices. Now, in Figure \ref{Fig}, we consider the Doily representation of a generalize quadrangle of order $(2,2)$ with nine marked flags (i.e., incident point-line pairs) with distinct non-black colours such that no two of them share a point or a line. Then the collection of these nine marked flags represents a perfect $1$-code in the line graph of the Tutte-Coxeter graph since each of its vertices corresponds to a flag of $\operatorname{GQ(2,2)}$. We note that the points of these nine marked flags together with the black lines form a $\operatorname{GQ(2,1)}$ subquadrangle, and their lines together with the black points form a complementary $\operatorname{GQ(1,2)}$ subquadrangle. Now we can conclude the following proposition.
\begin{center}
\begin{figure}[h]

 \centering

 \centerline{\includegraphics*[height=7cm]{Tutte-new}}
 \label{diagram}
 \caption{The Doily representation of $\operatorname{GQ(2,2)}$} \label{Fig}
\end{figure}
\end{center}
\begin{prop} \label{LGQ}
The line graph of the incidence graph of a generalized quadrangle of order $(q,q)$ admits a perfect $1$-code if and only if $q=2$.
\end{prop}
\subsection{The line graph of the incidence graph of a generalized hexagon} \label{SLGH}
Let $\Gamma$ be the incidence graph of a generalized hexagon of order $(q,q)$, where $q \geq 2$. Then it has intersection array $$\{q+1,q,q,q,q,q;1,1,1,1,1,q+1\}$$ and seven distinct eigenvalues $\{\pm(q+1),\pm \sqrt{3q},\pm \sqrt{q},0\}$. If the line graph of this graph admits a perfect $1$-code, then the only possible case is $q=3$ by Proposition \ref{bipartite}. Let $\Gamma$ be the incidence graph of a generalized hexagon of order $(3,3)$. If the line graph $\operatorname{L(\Gamma)}$ admits a perfect $1$-code $C$, then $|C|=208$ by Equation \ref{equation} since this graph is $6$-regular with $1456$ vertices. On the other hand, the code $C$ in the line graph $\operatorname{L(\Gamma)}$ can be related to a code $\overline{C}$ in the graph $\Gamma$ in such a way that $\overline{C}$ is a $1$-regular induced subgraph of the graph $\Gamma$ consisting of $208$ edges with mutually disjoint closed edge neighbourhoods. Additionally, $\overline{C}$ is also a vertex cover of the graph $\Gamma$.

Recall that the only known generalized hexagon of order $(3,3)$ is called the \emph{split Cayley hexagon of order $3$}. Let $\Gamma$ be the incidence graph of the split Cayley hexagon of order $3$. Then it has a lot of substructures which are generalized hexagons of order $(1,3)$ and $(3,1)$ (cf. \cite{DV}). Every generalized hexagon of order $(1,q)$ is isomorphic to the double of a projective plane of order $q$ for which the point set of the double is the set of points and lines and the line set is exactly the set of flags of the projective plane.
\begin{lem} \label{GH(1,q)}
The incidence graph of the double of a projective plane of order $q$ contains at most $q^{2}+q+1$ edges with mutually disjoint closed edge neighbourhoods.
\end{lem}
\begin{proof}
Let $\overline{C_{1}}$ be a collection of edges with mutually disjoint closed edge neighbourhoods of the incidence graph of this structure. If an edge $(p,(p,\ell))$ is in the set $\overline{C_{1}}$, for a point $p$ and line $\ell$, then $(\ell,(p',\ell))$ is not in $\overline{C_{1}}$ for every point $p'$ in $\ell$. Similarly, if an edge $(\ell,(p,\ell))$ is in the set $\overline{C_{1}}$, for a point $p$ and line $\ell$, then $(p,(p,\ell'))$ is not in $\overline{C_{1}}$ for every line $\ell'$ which contains the point $p$. On the other hand, a projective plane of order $q$ contains $q^{2}+q+1$ points and $q^{2}+q+1$ lines which implies that the set $\overline{C_{1}}$ contains at most $q^{2}+q+1$ edges with mutually disjoint closed edge neighbourhoods and the result follows.
\end{proof}

By using the {\sf FinInG} package \cite{fining} for {\sf GAP} \cite{GAP}, the graph $\Gamma$ can be constructed as the incidence graph of a block desig
\footnote{It has $364$ points and $364$ blocks such that each block consists of $4$ points and each point belongs to $4$ blocks.} with the point set $\{1,2,\ldots,364\}$ by the following commands.
\begin{verbatim}
gh:=SplitCayleyHexagon(3);
des:=BlockDesignOfGeneralisedPolygon(gh);
\end{verbatim}
Without loss of generality, let the point $1$ be outside the code $\overline{C}$. Then this vertex belongs to a substructure with the following point set $A$ of size $26$ and block set $B$ of size $52$ consisting of all blocks of this design with the property that each point in $A$ belongs to $4$ blocks in $B$.
\begin{verbatim}
A:={1,2,3,4,5,6,7,8,9,12,15,18,23,28,33,44,65,68,69,88,91,92,129,130,178,179};
B:={[1,2,111,112],[1,28,55,56],[1,65,127,128],[1,88,176,177],[2,3,224,227],
[2,4,225,228],[2,5,223,226],[3,23,45,46],[3,129,245,247],[3,178,311,313],
[4,18,34,35],[4,130,246,248],[4,179,292,312],[5,6,10,11],[5,7,13,14],
[5,8,16,17],[6,44,93,94],[6,91,189,194],[6,92,190,195],[7,33,70,71],
[7,68,140,145],[7,69,141,146],[8,9,21,22],[8,12,26,27],[8,15,31,32],
[9,88,186,191],[9,178,321,327],[9,179,281,293],[12,65,137,142],
[12,129,259,265],[12,130,260,266],[15,18,38,41],[15,23,49,52],
[15,28,59,62],[18,69,156,167],[18,92,203,214],[23,68,155,166],
[23,91,202,213],[28,33,78,83],[28,44,101,106],[33,130,222,279],
[33,163,178,331],[44,129,171,300],[44,179,208,338],[65,68,175,344],
[65,92,218,297],[68,179,241,355],[69,88,210,334],[69,129,244,357],
[88,91,161,306],[91,130,242,361],[92,178,239,356]}
\end{verbatim}
In this substructure, if we consider each block as a line, then it is isomorphic to the generalized hexagon of order $(1,3)
\footnote{We double-checked this with {\sf GAP} \cite{GAP}.}. Moreover, each block in $B$ consists of two points in $A$ and two points in the following set $P$ of size $104$.
\begin{verbatim}
P:={10,11,13,14,16,17,21,22,26,27,31,32,34,35,38,41,45,46,49,52,55,56,59,62,
70,71,78,83,93,94,101,106,111,112,127,128,137,140,141,142,145,146,155,156,
161,163,166,167,171,175,176,177,186,189,190,191,194,195,202,203,208,210,213,
214,218,222,223,224,225,226,227,228,239,241,242,244,245,246,247,248,259,260,
265,266,279,281,292,293,297,300,306,311,312,313,321,327,331,334,338,344,355,
356,357,361}
\end{verbatim}
Now consider another substructure with point set $M$ consisting of the points not in $A \cup P$ and block set $N$ consisting of the blocks not in $B$. In this new substructure, if we consider each block as a line, then it is indeed isomorphic to the interesting subgeometry with $234$ points and $312$ lines which is illustrated in \cite[Sec.~$3$]{DV}. By these structures and the properties of the code $\overline{C}$, we could obtain a contradiction as follows. Let $\Gamma_{i}(u)$ denote the set of vertices at distance $i$ from the vertex $u$ in the graph $\Gamma$, where $1\leq i \leq 6$. Recall that without loss of generality, the point $1$ is outside the code $\overline{C}$. We use \emph{inner} and \emph{outer} for elements in $\overline{C}$ and outside $\overline{C}$, respectively. Therefore the vertex $1$ is adjacent to $4$ inner blocks. Moreover, these $4$ blocks are adjacent to $4$ inner points and $8$ outer points in $\Gamma_{2}(1)$ forming the sets $\Gamma_{2}^{in}(1)$ and $\Gamma_{2}^{out}(1)$, respectively. We proceed this approach to find $\Gamma_{6}^{in}(1)$ and $\Gamma_{6}^{out}(1)$. The $4$ inner points of $\Gamma_{2}^{in}(1)$ are adjacent to $12$ outer blocks in $\Gamma_{3}(1)$ forming the set $\Gamma_{3}^{out}(1)$, and the $8$ outer points of $\Gamma_{2}^{out}(1)$ are adjacent to $24$ inner blocks in $\Gamma_{3}(1)$ forming the set $\Gamma_{3}^{in}(1)$. Moreover, the $12$ outer blocks of $\Gamma_{3}^{out}(1)$ are adjacent to $36$ inner points in $\Gamma_{4}(1)$ forming the set $\Gamma_{4}^{in}(1)$, and the $24$ inner blocks of $\Gamma_{3}^{in}(1)$ are adjacent to $24$ inner points and $48$ outer points in $\Gamma_{4}(1)$ forming the sets $\Gamma_{4}^{in'}(1)$ and $\Gamma_{4}^{out}(1)$, respectively. Furthermore, the $36$ inner points of $\Gamma_{4}^{in}(1)$ are adjacent to $36$ inner blocks and $72$ outer blocks in $\Gamma_{5}(1)$ forming the sets $\Gamma_{5}^{in}(1)$ and $\Gamma_{5}^{out}(1)$, respectively, and the $24$ inner points of $\Gamma_{4}^{in'}(1)$ are adjacent to $72$ outer blocks in $\Gamma_{5}(1)$ forming the set $\Gamma_{5}^{out'}(1)$, and the $48$ outer points of $\Gamma_{4}^{out}(1)$ are adjacent to $144$ inner blocks in $\Gamma_{5}(1)$ forming the set $\Gamma_{5}^{in'}(1)$. Finally, the $144$ inner blocks of $\Gamma_{5}^{in'}(1)$ are adjacent to $144$ inner points in $\Gamma_{6}(1)$ forming the set $\Gamma_{6}^{in}(1)$ and the $99$ remaining points in $\Gamma_{6}(1)$ are outer points forming the set $\Gamma_{6}^{out}(1)$ since the code $\overline{C}$ consists of $4+24+36+144=208$ points.

By using {\sf GAP} \cite{GAP}, it turns out that there are $4$ points in $A \cap \Gamma_{2}(1)$ and $8$ points in $P \cap \Gamma_{2}(1)$. Moreover, there are $12$ blocks in $B \cap \Gamma_{3}(1)$ and $24$ blocks in $N \cap \Gamma_{3}(1)$. Furthermore, there are $12$ points in $A \cap \Gamma_{4}(1)$, $24$ points in $P \cap \Gamma_{4}(1)$ and $72$ points in $M \cap \Gamma_{4}(1)$. Moreover, there are $36$ blocks in $B \cap \Gamma_{5}(1)$ and $288$ blocks in $N \cap \Gamma_{5}(1)$. Finally, there are $9$ points in $A \cap \Gamma_{6}(1)$, $72$ points in $P \cap \Gamma_{6}(1)$ and $162$ points in $M \cap \Gamma_{6}(1)$.

Moreover, there exist five cases based on the number of inner points $i$ in $A \cap \Gamma_{2}(1)$, where $0 \leq i \leq 4$. If $|A \cap \Gamma_{2}^{in}(1)|=i$, then $|A \cap \Gamma_{2}^{out}(1)|=|P \cap \Gamma_{2}^{in}(1)|=4-i$, $|P \cap \Gamma_{2}^{out}(1)|=4+i$, $|B \cap \Gamma_{3}^{in}(1)|=|A \cap (\Gamma_{4}^{in'}(1) \cup \Gamma_{4}^{out}(1))|=|N \cap \Gamma_{3}^{out}(1)|=12-3i$, $|N \cap \Gamma_{3}^{in}(1)|=|M \cap \Gamma_{4}^{in'}(1)|=12+3i$, $|B \cap \Gamma_{3}^{out}(1)|=|A \cap \Gamma_{4}^{in}(1)|=|B \cap \Gamma_{5}^{in}(1)|=3i$, $|P \cap \Gamma_{4}^{in}(1)|=|B \cap \Gamma_{5}^{out}(1)|=6i$, $|P \cap (\Gamma_{4}^{in'}(1) \cup \Gamma_{4}^{out}(1))|=24-6i$, $|M \cap \Gamma_{4}^{in}(1)|=36-9i$, and $|M \cap \Gamma_{4}^{out}(1)|=24+6i$. Moreover, as each point of $P \cap \Gamma_{4}(1)$ has the neighbor from $B$ in $\Gamma_{3}(1)$, it follows that each block of $B \cap \Gamma_{5}(1)$ has one neighbor in $A \cap \Gamma_{4}(1)$ and one in $A \cap \Gamma_{6}(1)$.

Now suppose that $i \geq 1$. Then $|A \cap \Gamma_{6}^{in}(1)|=6$ and $|A \cap \Gamma_{6}^{out}(1)|=3$. Since each point of $A \cap \Gamma_{6}^{out}(1)$ has precisely $i$ neighbors in $B \cap \Gamma_{5}^{in}(1)$, it follows that $|B \cap \Gamma_{5}^{in'}(1)|=6+3(4-i)=18-3i$, and then $|A \cap \Gamma_{4}^{out}(1)|=6-i$ and $|A \cap \Gamma_{4}^{in'}(1)|=6-2i$. Therefore, there are $i+3i+(6-2i)+6=12+2i>13$ inner points in $A$, each of which is adjacent to an inner block in $B$, contradicting Lemma \ref{GH(1,q)} (see Figure \ref{Fig2}). This then leaves us with the case $i=0$.

Let there exist $a$ inner and $9-a$ outer points in $A \cap \Gamma_{6}(1)$, where $0 \leq a \leq 9 $. Then there are $3a$ blocks in $B \cap \Gamma_{5}^{out'}(1)$. Moreover, $|A \cap \Gamma_{4}^{in'}(1)|=a$ since each block of $B \cap \Gamma_{5}^{out'}(1)$ has one neighbor in $A \cap \Gamma_{4}^{in'}(1)$, and each point of $A \cap \Gamma_{4}^{in'}(1)$ has three neighbors in $B \cap \Gamma_{5}^{out'}(1)$. If $a \geq 1$, then consider the $(3,2)$-biregular bipartite incidence graph with the point set consisting of the $2a$ inner points in the union of $A \cap \Gamma_{4}^{in'}(1)$ and $A \cap \Gamma_{6}(1)$, and the block set consisting of the $3a$ blocks in $B \cap \Gamma_{5}^{out'}(1)$. This graph contains at least $35$ vertices since the girth of this graph is at least $12$. This implies that $a \geq 7$ and therefore there are at least $14$ inner points in $A$, contradicting Lemma \ref{GH(1,q)}. It follows that the only possible case for $i=0$ is as in Figure \ref{Fig3} and we prove that it is impossible. On one side, each of the $12$ points in $P \cap \Gamma_{4}^{in'}(1)$ has three neighbors in $N \cap \Gamma_{5}^{out'}(1)$ and two points in the set $P$ can not share the same neighbor in $N$. Therefore each block in the set $R$, consisting of the $36$ blocks in $N \cap \Gamma_{5}^{out'}(1)$ which do not have a neighbor in $P \cap \Gamma_{4}^{in'}(1)$, must be adjacent to a unique point in $P \cap \Gamma_{6}^{in}(1)$. Moreover, each of the $36$ points in $P \cap \Gamma_{6}^{in}(1)$ has exactly one neighbor in $N \cap \Gamma_{5}^{out'}(1)$.  To see this, suppose in contrary that there exists a point in $P \cap \Gamma_{6}^{in}(1)$ which has at least two neighbors in $N \cap \Gamma_{5}^{out'}(1)$. Then there exists a point $u \in P \cap \Gamma_{6}^{in}(1)$ which has no neighbor in $R$. Furthermore, there are two blocks in $B \cap \Gamma_{5}^{in'}(1)$ at distance $3$  and therefore two points in $P \cap \Gamma_{6}^{in}(1)$ at distance $4$ form $u$ which are adjacent to at most six blocks in $R$. Moreover, there are at most three blocks in $N \cap \Gamma_{5}^{out}(1)$ at distance $3$ and therefore at most three points in $P \cap \Gamma_{6}^{in}(1)$ at distance $4$ form $u$ which are adjacent to at most six blocks in $R$. Additionally, there are no blocks in $N \cap \Gamma_{3}^{in}(1)$ at distance $3$ from $u$. It follows that there is no path of length at most $6$ from $u$ to some blocks of $R$, a contradiction. This implies that each of the $36$ points in $P \cap \Gamma_{6}^{in}(1)$ has exactly two neighbors in $N \cap \Gamma_{5}^{out}(1)$. On the other side, by using {\sf GAP} \cite{GAP}, it turns out that there are exactly $27$ points in $P \cap \Gamma_{6}^{in}(1)$ which have two neighbors in $\Gamma_{5}^{in}(1) \cup \Gamma_{5}^{out}(1)
\footnote{There are $2^{4}=16$ cases for these sets depending on the choice of the four inner points in $P \cap \Gamma_{2}^{in}(1)$. All of these cases have been checked with {\sf GAP} \cite{GAP}.}, a contradiction, and this completes the proof.

\begin{center}
\begin{figure}[h]

 \centering

 \centerline{\includegraphics*[height=12cm]{GHStructure2}}
 \label{diagram}
 \caption{The incidence graph of $\operatorname{GH(3,3)}$ for $i>0$} \label{Fig2}
\end{figure}
\end{center}

\begin{center}
\begin{figure}[h]

 \centering

 \centerline{\includegraphics*[height=12cm]{GHStructure}}
 \label{diagram}
 \caption{The incidence graph of $\operatorname{GH(3,3)}$ for $i=0$} \label{Fig3}
\end{figure}
\end{center}

Therefore we can conclude the following proposition.
\begin{prop} \label{LGH}
If the line graph of the incidence graph of a generalized hexagon of order $(q,q)$, where $q \geq 2$, admits a  perfect $1$-code, then $q=3$. Moreover, the incidence graph the split Cayley hexagon of order $3$ doesn't admit a perfect $1$-code
\footnote{The anonymous referee double-checked this with {\sf GAP} \cite{GAP}.}
\end{prop}
\section{Distance-regular graphs with small diameter}
As far as we know, the general problem of characterizing distance-regular graphs with small diameter greater than $2$ which admit a perfect $1$-code is hard. Therefore we give an overview up to diameter $4$. For a complete graph, the only perfect $1$-codes are the one-element subsets (cf. Prop.~\ref{completegraph}). Moreover, by Observation \ref{Observation1}, there is no perfect $1$-code in a strongly regular graph.

Now suppose that $C$ is a perfect $1$-code in a distance-regular graph with diameter $3$. Then the distance between two vertices in $C$ must be $3$.
Let $C$ be a perfect $1$-code in a bipartite distance-regular graph with diameter $3$. Then $C$ contains exactly two elements from different parts since the distance between two vertices in $C$ must be $3$. Therefore this graph must be a complete bipartite graph minus a perfect matching. Hence we can conclude the following proposition.
\begin{prop} \label{Bipartitedrg3}
Let $\Gamma$ be a bipartite distance-regular graph with diameter $3$. Then $\Gamma$ admits a perfect $1$-code $C$ if and only if $\Gamma$ is a complete bipartite graph minus a perfect matching and $C$ contains exactly two elements at distance $3$ from different parts.
\end{prop}
Note that every antipodal distance-regular graph of diameter $3$ (including the bipartite ones) admits a perfect $1$-code -- in fact, each fiber is such a code. Now we deal with primitive distance-regular graphs with diameter $3$.
\begin{Observation} \label{Observation3}
If a primitive distance-regular graph $\Gamma$ with diameter $3$ admits a perfect $1$-code, then it has eigenvalue $-1$ by Observation \ref{Observation2} and therefore its distance-$3$ graph is strongly regular (cf. \cite[Prop.~4.2.17]{BCN}). It follows that if this strongly regular graph has parameters $(n,k,\lambda,\mu)$, then $\lambda \geq |C|-2$ because the perfect code $C$ is a clique in the distance-$3$ graph by Observation \ref{Observation1}.
\end{Observation}
Among primitive distance-regular graphs with diameter $3$ and small number of vertices, the first putative example is the Odd graph with $35$ vertices and intersection array $\{4,3,3;1,1,2\}$ since it has eigenvalue $-1$ (cf. \cite[Chap.~$14$]{BCN} and Obs. \ref{Observation2}). This graph admits a perfect $1$-code with $7$ vertices (cf. \cite[Fig.~1]{JV}). The second example is the Sylvester graph with $36$ vertices and intersection array $\{5,4,2;1,1,4\}$. Indeed every $6$-clique in the distance-$3$ graph of the Sylvester graph corresponds to a perfect $1$-code in this graph (see also \cite[Sec.~$3$]{JV2}).
\begin{prop} \label{Sylvester}
The Sylvester graph admits a perfect $1$-code.
\end{prop}

Let $\Gamma$ be a distance-regular graph with diameter $4$. If $\Gamma$ is antipodal, then it doesn't admit a perfect $1$-code (cf. Cor. \ref{antipodaldrg45}). Now suppose that $\Gamma$ is bipartite with degree $k$. Then it has eigenvalues $\{\pm k,0,\pm1\}$ because it has eigenvalue $-1$ by Observation \ref{Observation2}. The following lemma shows there is no such graph.
\begin{lem}
There is no bipartite distance-regular graph with diameter $4$ and eigenvalues $\{\pm k,0,\pm1\}$.
\end{lem}
\begin{proof}
Let $\Gamma$ be a bipartite distance-regular graph with diameter $4$ and eigenvalues $\{\pm k,0,\pm1\}$. We note that this graph can not be the cycle graph on $8$ vertices and therefore $k>2$. On the other hand, if this graph has intersection array $\{k,b_{1},b_{2},b_{3};1,c_{2},c_{3},c_{4}\}$, then by considering the intersection matrix of this graph, the second largest eigenvalue must be $(c_{2}+1)k-c_{2}(c_{3}+1)$ and therefore $(c_{2}+1)k-c_{2}(c_{3}+1)=1$. On the other hand, $c_{2}$ must divide $k$ (cf. \cite[Lem.~$1.7.2$]{BCN}) which implies that $c_{2}=1$ and $c_{3}=2k-2$. Therefore $k \leq 2$ since $c_{3} \leq k$, a contradiction, and this completes the proof.
\end{proof}
Therefore we can conclude the following proposition.
\begin{prop}
A bipartite distance-regular graph with diameter $4$ doesn't admit a perfect $1$-code.
\end{prop}
Among primitive distance-regular graphs with diameter $4$ and small number of vertices, the first putative example is the Coxeter graph with $28$ vertices and intersection array $\{3,2,2,1;1,1,1,2\}$ since it has eigenvalue $-1$ (cf. \cite[Chap.~$14$]{BCN} and Obs. \ref{Observation2}).
This graph partitions into three $7$-gons and a $7$-coclique (cf. \cite{Coxeter}). It turns out that the $7$-coclique is indeed a perfect $1$-code. Therefore we can conclude the following proposition.
\begin{prop} \label{Coxeter}
The Coxeter graph admits a perfect $1$-code.
\end{prop}
\section{Distance-regular graphs with small valency}
All known distance-regular graphs with small valency at most $4$, the distance-regular graphs with known putative intersection arrays for valency $5$, and all distance-regular graphs with girth $3$ and valency $6$ or $7$  are listed in \cite{VJ}. We give an overview of all possible intersection arrays and corresponding graphs, and indicate which of these admit a perfect $1$-code. Note that for each intersection array in Table \ref{tabledrgvalency3} there is a unique distance-regular graph. Moreover, for each intersection array in Table \ref{tabledrgvalency4} there is a unique distance-regular graph, except possibly for the last array, which corresponds to the incidence graph of a generalized hexagon of order $(3,3)$. Furthermore, in Table \ref{tabledrgvalency5}, all known putative intersection arrays for distance-regular graphs with valency $5$ are listed.  All of the graphs in the table are unique, given their intersection arrays, except possibly the incidence graph of a generalized hexagon of order $(4,4)$ (the last case). All distance-regular graphs with valency at most $7$ and girth $3$ (i.e., with triangles) are listed in Table \ref{tabledrgvalency6} besides the ones with valency at most $5$ that we have encountered in the previous tables. For each of the intersection arrays $\{6,3;1,2\}$ and $\{6,4,4;1,1,3\}$, there are exactly two distance-regular graphs (as mentioned in the table). By $n$, $d$, and $g$, we denote the number of vertices, diameter, and girth, respectively. We note that in the reference column, only one reason is stated.
\subsection{The point graphs of the generalized hexagons of order $(2,2)$} \label{GH(2,2)}
Up to isomorphism there are exactly two generalized hexagons of order $(2,2)$. Each of them is the dual of the other (cf. \cite[Theorem~1]{CT}). Their point graphs (collinearity graphs) give rise to two distance-regular graphs with intersection array $\{6,4,4;1,1,3\}$ and $63$ vertices. We can distinguish the two graphs by whether the graph induced on the vertices at distance $3$ from a fixed vertex is connected or not. These two non-isomorphic distance-regular graphs have been constructed in {\sf GAP} \cite{GAP} with {\sf Grape} \cite{Sgrape} package as {\sf Graph $1$} and {\sf Graph $2$} in \cite{drgorg}. Indeed the graph induced on the vertices at distance $3$ from a fixed vertex in the {\sf Graph $2$} is connected whereas in the {\sf Graph $1$} is disconnected. Let $\Gamma$ be one of such distance-regular graphs which admits a perfect $1$-code $C$. Then $|C|=9$ by Equation \ref{equation}. Moreover, each pair of vertices in $C$ is at distance $3$ by Observation \ref{Observation1} since $\Gamma$ has diameter $3$. Therefore a perfect $1$-code $C$ can be viewed as a clique in the distance-$3$ graph of $\Gamma$. By using {\sf GAP} \cite{GAP}, it turns out that the {\sf Graph $1$} of \cite{drgorg} admits perfect $1$-code
\footnote{The {\sf Graph $1$} of \cite{drgorg} with vertex set $\{1,2,\ldots,63\}$ admits $\{3,4,5,7,35,37,42,50,63\}$ and $\{1,2,22,24,30,33,41,61,63\}$ as perfect $1$-codes.} but the {\sf Graph $2$} of \cite{drgorg} doesn't admit a perfect $1$-code since there is no complete subgraph of size $9$ in its distance-$3$ graph. This shows that we can not deduce that a distance-regular graph admits a perfect $1$-code from its spectrum.

\begin{center}
\begin{table}[h]
\begin{tabular}{ l r r c l  c l}
  \hline
  Intersection array & $n$ & $d$ & $g$ & Name & Perfect $1$-code & Reference\\
  \hline

  \{3;1\}& 4 & 1 & 3&$\operatorname{K_{4}}$ & Yes & Rem.~\ref{completegraph}\\
  \{3,2;1,3\} & 6 & 2 &4& $\operatorname{K_{3,3}}$ & No & Eq.~\ref{equation}\\
  \{3,2,1;1,2,3\} & 8 & 3 &4& $\operatorname{K_{3,3}^*}$ & Yes & Prop.~\ref{Bipartitedrg3}\\
  \{3,2;1,1\} & 10 & 2 &5& $\operatorname{Petersen}$ & No & Eq.~\ref{equation}\\
  \{3,2,2;1,1,3\} & 14 & 3 &6& $\operatorname{Heawood}$ & No & Eq.~\ref{equation}\\
  \{3,2,2,1;1,1,2,3\} & 18 & 4 &6& $\operatorname{Pappus}$ & No & Eq.~\ref{equation}\\
  \{3,2,2,1,1;1,1,2,2,3\} & 20 & 5 &6& $\operatorname{Desargues}$ & No & Cor.~\ref{antipodaldrg45}\\
  \{3,2,1,1,1;1,1,1,2,3\} & 20 & 5 &5& $\operatorname{Dodecahedron}$ & No & Cor.~\ref{antipodaldrg45}\\
  \{3,2,2,1;1,1,1,2\}& 28 & 4 &7& $\operatorname{Coxeter}$ & Yes & Prop.~\ref{Coxeter} \\
  \{3,2,2,2;1,1,1,3\} & 30 & 4 &8& $\operatorname{Tutte's}$ $\operatorname{8-cage}$ & No & Eq.~\ref{equation}\\
  \{3,2,2,2,2,1,1,1;& 90 & 8 & 10& $\operatorname{Foster}$ & No & Eq.~\ref{equation}\\
   \hspace{.5cm}1,1,1,1,2,2,2,3\}  & && && &\\
  \{3,2,2,2,1,1,1; & 102 & 7 &9& $\operatorname{Biggs-Smith}$ & No & Eq.~\ref{equation}\\
  \hspace{.5cm}1,1,1,1,1,1,3\}  && & & &&\\
  \{3,2,2,2,2,2; & 126 & 6 &12& $\operatorname{Tutte's}$ $\operatorname{12-cage}$ & No & Eq.~\ref{equation}\\
   \hspace{.5cm}1,1,1,1,1,3\}  & && & &&\\
  \hline
\end{tabular}
\caption{Distance-regular graphs with valency $3$}\label{tabledrgvalency3}
\end{table}
\end{center}
\begin{table}[h]
\begin{center}
\begin{tabular}{ l r r c l c l}
\hline
   Intersection array & $n$ & $d$ & $g$& Name & Perfect $1$-code & Reference\\
  \hline

  \{4;1\} & 5 & 1 &3& $\operatorname{K_{5}}$ & Yes & Rem.~\ref{completegraph}\\
  \{4,1;1,4\} &  6 & 2 &3& $\operatorname{K_{2,2,2}}$ & No & Eq.~\ref{equation}\\
  \{4,3;1,4\} &  8 & 2 &4& $\operatorname{K_{4,4}}$ & No & Eq.~\ref{equation}\\
  \{4,2;1,2\} &  9 & 2 &3& $\operatorname{Paley}$ $\operatorname{graph}$ $\operatorname{P(9)}$ & No & Eq.~\ref{equation}\\
  \{4,3,1;1,3,4\} &  10 & 3 &4& $\operatorname{K^{*}_{5,5}}$ & Yes & Prop.~\ref{Bipartitedrg3}\\
  \{4,3,2;1,2,4\} &  14 & 3 &4& $\operatorname{IG(7,4,2)}$ & No & Eq.~\ref{equation}\\
  \{4,2,1;1,1,4\} &  15 & 3 &3& $\operatorname{L(Petersen)}$ & Yes & Cor.~\ref{antipodaldrg3}\\
  \{4,3,2,1;1,2,3,4\} &  16 & 4 &4& $\operatorname{Q_{4}}$ & No & Eq.~\ref{equation}\\
  \{4,2,2;1,1,2\} &  21 & 3 &3& $\operatorname{L(Heawood)}$ & No & Eq.~\ref{equation}\\
   \{4,3,3;1,1,4\} &  26 & 3 &6& $\operatorname{IG(13,4,1)}$ & No & Eq.~\ref{equation}\\
   \{4,3,3,1;1,1,3,4\} &  32 & 4 &6& $\operatorname{IG(A(2,4)\setminus pc)}$ & No & Eq.~\ref{equation}\\
   \{4,3,3;1,1,2\} &  35 & 3 &6& $\operatorname{O_4}$ & Yes & \cite[Fig.~1]{JV}\\
   \{4,2,2,2;1,1,1,2\} &  45 & 4 &3& $\operatorname{L(Tutte's}$ $\operatorname{8-cage)}$
   & Yes & Prop.~\ref{LGQ}\\%\cite[Prop.~3.4]{ADJ}\\
   \{4,3,3,2,2,1,1; &  70 & 7 &6& $\operatorname{DO_4}$ & Yes & Cor.~\ref{DoubleOdd}\\
    \hspace{.5cm}1,1,2,2,3,3,4\} & & & && &\\
   \{4,3,3,3;1,1,1,4\} &  80 & 4 &8& $\operatorname{IG(GQ(3,3))}$ & No & Obs.~\ref{Observation2}\\
   \{4,2,2,2,2,2; &  189 & 6 &3& $\operatorname{L(Tutte's}$ $\operatorname{12-cage})$
   & No & Eq.~\ref{equation}\\
     \hspace{.5cm}1,1,1,1,1,2\} & & & && &\\
   \{4,3,3,3,3,3; &  728 & 6 &12& $\operatorname{IG(GH(3,3))}$ & No & Eq.~\ref{equation}\\
     \hspace{.5cm}1,1,1,1,1,4\}  && & && &\\
  \hline
\end{tabular}
 \caption{Distance-regular graphs with valency $4$}\label{tabledrgvalency4}
\end{center}
\end{table}
\begin{table}[h]
\begin{center}
\begin{tabular}{ l r r c l c l}
  \hline
   Intersection array & $n$ & $d$ &$g$& Name & Perfect $1$-code & Reference\\
  \hline

 \{5;1\} & 6 & 1 &3& $\operatorname{K_{6}}$ & Yes & Rem.~\ref{completegraph}\\
 \{5,4;1,5\} &  10 & 2 &4& $\operatorname{K_{5,5}}$ & No & Eq.~\ref{equation}\\
 \{5,2,1;1,2,5\} &  12 & 3 &3& $\operatorname{Icosahedron}$ & Yes & Cor.~\ref{antipodaldrg3}\\
 \{5,4,1;1,4,5\} &  12 & 3 &4 &$\operatorname{K^{*}_{6,6}}$ & Yes & Prop.~\ref{Bipartitedrg3}\\
 \{5,4;1,2\} &  16 & 2 &4& $\operatorname{Folded}$ $\operatorname{5-cube}$ & No & Eq.~\ref{equation}\\
 \{5,4,3;1,2,5\} &  22 & 3 &4&$\operatorname{IG(11,5,2)}$ & No & Eq.~\ref{equation}\\
 \{5,4,3,2,1;1,2,3,4,5\} &  32 & 5 &4& $\operatorname{Q_{5}}$ & No & Eq.~\ref{equation}\\
 \{5,4,1,1;1,1,4,5\} &  32 & 4 &5& $\operatorname{Armanios-Wells}$ & No & Eq.~\ref{equation}\\
 \{5,4,2;1,1,4\} &  36 & 3 &5& $\operatorname{Sylvester}$ & Yes & Prop.~\ref{Sylvester}\\
  \{5,4,4;1,1,5\} &  42 & 3 &6& $\operatorname{IG(21,5,1)}$ & No & Obs.~\ref{Observation2}\\
  \{5,4,4,1;1,1,4,5\} &  50 & 4 &6& $\operatorname{IG(A(2,5)\setminus pc)}$ & No & Eq.~\ref{equation}\\
  \{5,4,4,3;1,1,2,2\} &  126 & 4 &6& $\operatorname{O_{5}}$ & No & Obs.~\ref{Observation2}\\
  \{5,4,4,4;1,1,1,5\} &  170 & 4 &8& $\operatorname{IG(GQ(4,4))}$ & No & Eq.~\ref{equation}\\
  \{5,4,4,3,3,2,2,1,1; &  252 & 9 &6& $\operatorname{DO_{5}}$ & No & Cor.~\ref{DoubleOdd}\\
   \hspace{.5cm}1,1,2,2,3,3,4,4,5\}  && & && &\\
  \{5,4,4,4,4,4; &  2730 & 6 &12& $\operatorname{IG(GH(4,4))}$ & No & Obs.~\ref{Observation2}\\
  \hspace{.5cm}1,1,1,1,1,5\}  && & && &\\
  \hline
\end{tabular}
 \caption{Distance-regular graphs with valency $5$}\label{tabledrgvalency5}
\end{center}
\end{table}
\begin{table}[h]
\begin{center}
\begin{tabular}{ l r r c l c l}
  \hline
   Intersection array & $n$ & $d$ &$g$& Name & Perfect $1$-code & Reference\\
  \hline

 \{6;1\} & 7 & 1 &3& $\operatorname{K_{7}}$ & Yes & Rem.~\ref{completegraph}\\
 \{6,1;1,6\} &  8 & 2 &3& $\operatorname{K_{2,2,2,2}}$ & No & Eq.~\ref{equation}\\
 \{6,2;1,6\} &  9 & 2 &3& $\operatorname{K_{3,3,3}}$ & No & Eq.~\ref{equation}\\
 \{6,2;1,4\} &  10 & 2 &3& $\operatorname{T(5)}$ & No & Eq.~\ref{equation}\\
 \{6,3;1,3\} &  13 & 2 &3& $\operatorname{P(13)}$ & No & Eq.~\ref{equation}\\
 \{6,4;1,3\} &  15 & 2 &3& $\operatorname{\overline{T(6)}\sim GQ(2,2)}$ & No & Eq.~\ref{equation}\\
 \{6,3;1,2\} &  16 & 2 &3& $\operatorname{L_{2}(4)}$, $\operatorname{Shrikhande}$ & No & Eq.~\ref{equation}\\
 \{6,4,2;1,2,3\} &  27 & 3 &3& $\operatorname{H(3,3)}$ & No & Eq.~\ref{equation}\\
 \{6,4,2,1;1,1,4,6\} &  45 & 4 &3& $\operatorname{halved}$ $\operatorname{Foster}$ & No & Eq.~\ref{equation}\\
 \{6,3,3;1,1,2\} &  52 & 3 &3& $\operatorname{L(IG(13,4,1)}$ & No & Eq.~\ref{equation}\\
 \{6,4,4;1,1,3\} &  63 & 4 &3& $\operatorname{GH(2,2)}$ $\operatorname{(Graph~1)}$& Yes & Sec.~\ref{GH(2,2)}\\
 \{6,4,4;1,1,3\} &  63 & 4 &3& $\operatorname{GH(2,2)}$ $\operatorname{(Graph~2)}$& No & Sec.~\ref{GH(2,2)}\\
 \{6,3,3,3;1,1,1,2\} &  160 & 4 &3& $\operatorname{L(IG(GQ(3,3)))}$ & No & Prop.~\ref{LGQ}\\
 \{6,3,3,3,3,3;1,1,1,1,1,2\} &  1456 & 6 &3& $\operatorname{L(IG(GH(3,3)))}$ & No & Prop.~\ref{LGH}\\
 \hline
  \{7;1\} & 8 & 1 &3& $\operatorname{K_{8}}$ & Yes & Rem.~\ref{completegraph}\\
 \{7,4,1;1,2,7\} &  24 & 3 &3& $\operatorname{Klein}$ & Yes & Cor.~\ref{antipodaldrg3}\\
 \hline
\end{tabular}
 \caption{Distance-regular graphs with girth $3$ and valency $6$ or $7$}\label{tabledrgvalency6}
\end{center}
\end{table}
\section*{Acknowledgements}

\noindent The author would like to thank the anonymous referee for his/her invaluable comments which led to fixing some errors in the proof of Proposition \ref{LGH} and the statement of Proposition \ref{Sylvester}, and improving the presentation of this paper. The author is grateful to the Research Council of Shahid Chamran University of Ahvaz for financial support (SCU.MM99.29248).


\section{\textsc{Introduction}}

The rise of hyperbolic materials in recent years promises important applications in
optoelectronics and nanophotonics~\cite{hyperbolic1,hyperbolic2,hyperbolic3,hyperbolic4,CuS19}. Light can acquire
hyperbolic dispersion while passing through such materials, which occurs in some frequency range, when different principal components of the longitudinal dielectric function (dielectric permittivity) have opposite signs. For the case of isotropic medium, it behaves as a dielectric, that is, supports propagating electromagnetic waves, when the sign of the dielectric function is positive. When the latter is negative, the incident light is reflected, with only an exponentially decaying evanescent field penetrating the material, like for metals below plasma threshold.

Anisotropy in electronic, optical, vibrational, and transport behaviour can occur when structural anisotropy is present, and if it is sufficiently strong, the different components of the dielectric permittivity tensor may aquire opposite signs to turn the material hyperbolic.  In two-dimensional crystals, in-plane anisotropy strong enough to make it hyperbolic is a unique situation that allows to confine short wavelengths (large wave vectors) inside a material, promising smaller sizes for optoelectronic devices.  Optical anisotropy in rhenium disulfide ($\text{ReS}_2$) has been established both for bulk crystals \cite{lin2011anisotropy,friemelt1993optical} and for thin layers \cite{liu2015integrated}.  In this work we predict that ($\text{ReS}_2$) that appears in a distorted 1T phase can realize hyperbolic plasmons depending on the number of layers.

It has been suggested that anisotropic 2D materials can be tuned to become hyperbolic via electrostatic tuning, strain or dimensionality and can host hyperbolic plasmons \cite{nemilentsau2016anisotropic,van2019tuning}. Several studies \cite{low2017polaritons,sun2014indefinite,cordova2019anisotropic} have investigated hyperbolic plasmons (HP) and its existence in naturally occurring materials. The strong anisotropy of $\text{ReS}_2$ hints at its potential as a natural hyperbolic material, offering possibilities for studying HP. Previous works
\cite{echeverry2018theoretical,zhang2015res2} have studied the band structure and anisotropic optical response of
$\text{ReS}_2$, but the study of HP remains largely unexplored. Here we use the ladder-vertex corrected and local-field corrected plasmonic response in ReS$_{2}$ within a self-consistent solution of Bethe-Salpeter equation (BSE) as implemented in \questaal \cite{pashov2020questaal}.

The rest of the paper is organized as follows. In Sec. \ref{sec: atomic-str}, we describe the anisotropic atomic structure of $\text{ReS}_2$ and in Sec. \ref{sec: comp-details}, we briefly describe theoretical methods and provide computational details. In Secs. \ref{sec: band-str} and \ref{sec: optical-spectra} we present our
results on the electronic structure and optical properties of $\text{ReS}_2$. In Sec. \ref{sec: conclusion}, we briefly summarize our results and conclude the paper.

\section{\textsc{Atomic Structure and Computational Details}}
\subsection{\textsc{Atomic structure}} \label{sec: atomic-str}

\begin{figure*}[htb]
	\includegraphics[scale=0.3]{res2-str.pdf}%
	\caption{The ball-stick model of the bulk distorted 1T diamond-chain $\text{ReS}_2$ obtained using VESTA software \cite{momma2011vesta} is presented. Top view (Fig.1(a)) and side view (Fig.1(b))of the distorted 1T-$\text{ReS}_2$; Re atoms are in red and S atoms are in green. The black outline shows the unit cell used for the calculation. The Re chain is along b direction. Brillouin zone (bottom) of the corresponding hexagonal lattice with lines connecting high-symmetry points $\Gamma$-K1-K2-M2-$\Gamma$-K4-K3-M2-$\Gamma$. $a^*$ and $b^*$ denote reciprocal
		lattice vectors.}
	\label{fig:res2-str}

\end{figure*}

$\text{ReS}_2$ belongs to the family of two-dimensional (2D) layered transition metal dichalcogenide (TMDs) of the form $\text{MX}_2$ where M is a transition metal atom (Mo, W, Re, ...) and X is a group-16 atom (S, Se, Te). The atomic structure of $\text{ReS}_2$ layers has neither H or T character. Unlike other TMDs, which usually have 1H or 1T structure in their ground state, $\text{ReS}_2$ crystallizes in a distorted-1T structure with clustering of Re units forming parallel metal chains along the van der Waals plane (see Fig. \ref{fig:res2-str}).

The compound $\text{ReS}_2$ belongs to the triclinic symmetry group $\text{P}^1$, resembling a distorted CdCl$_{2}$ structure. It comprises three atomic layers, S-Re-S, where covalent bonds join Re and S. The adjacent layers of $\text{ReS}_2$ are coupled by weak van der Waals (vdW) forces to form bulk crystals. The unit cell is derived from hexagonal symmetry towards a distorted 1T structure, in which Re atoms group into parallelograms of four Re atoms. The formation of Re chains breaks the hexagonal symmetry and doubles the unit cell size.  Hence the unit cell of single layer $\text{ReS}_2$ in the distorted-1T phase is composed of four Re and eight S atoms. In pristine $\text{ReS}_2$, the valence band maximum is composed from 5d orbitals of Re atoms and 3p orbitals of S atoms, and the conduction band minimum is derived from 5d orbitals of Re atoms. The Brillouin zone of the $\text{ReS}_2$ is hexagonal but with unequal sides as a result of the distorted atomic structure.  Energy band structures were generated along the symmetry lines shown in Fig. \ref{fig:res2-str}.

\subsection{\textsc{Computational Details} } \label{sec: Computational_details}\label{sec: comp-details}

\subsubsection{\textbf{LDA, QS\emph{GW}, and \bse self-consistency}}

Single-particle calculations (LDA, and the the quasi-particle self-consistent \emph{GW}~\cite{van2006quasiparticle}
(QS\emph{GW}) self-energy $\Sigma^0(k)$) were performed on a $12 \times 12 \times 12$ points (Monkhorst pack) for bulk
and $12 \times 12 \times 1$ for ML and BL. An energy cutoff of 400 \si{\electronvolt} was used and Gaussian smearing
with a width of 0.05 \si{\electronvolt}. The tolerance of $10^{-5}$ \si{\electronvolt} and $2 \times 10^{-6}$ has been
taken in convergence of energy and RMS density respectively. The charge density was made self-consistent for each
iteration in the QS\emph{GW} self-consistency cycle. The QS\emph{GW} cycle was iterated until the RMS change in
$\Sigma^0$ reached $10^{-5}$ Ry. Thus the calculation was self-consistent in both $\Sigma^0(k)$ and the
density. Numerous checks were made to verify that the self-consistent $\Sigma^0(k)$(k) was independent of starting
point.  For ML-$\text{ReS}_2$, we performed a rigorous check for vacuum correction to all band gap and dielectric
screening by increasing the size from \SI{10}{\angstrom} to \SI{45}{\angstrom}. Since along the \emph{z}-direction we
have a vacuum, the dielectric constant, which is the real part of the macroscopic dielectric response at $\omega{=}0$,
should be close to unity.

In the present work, the electron-hole two-particle correlations are incorporated within a
self-consistent ladder BSE implementation \cite{cunningham2018effect} with Tamn-Dancoff approximation
\cite{hirata1999time}. Ladder diagrams are included in the polarizability $P$ that makes $W$, via the
  solution of a Bethe-Salpeter equation (BSE); thus this form of $GW$ goes beyond the RPA in constructing the
  self-energy $\Sigma=iGW$. The electron-hole attraction from the ladders enhances $P$, thus reducing $W$, which in turn
  reduces the bandgap.  A static vertex is used to construct $P$.  \emph{G} and \emph{W} are calculated
  self-consistently, in quasiparticlized form~\cite{van2006quasiparticle}: \emph{G} and \emph{W} are updated iteratively
  until all of them converge (QS\emph{GW}).  When ladders are incorporated into \emph{W}, we denote the process as
  QS$G\widehat{W}$ to signify $W$ was computed from the BSE.  The macroscopic dielectric function we present here,
  $[\epsilon^{-1}_{\mathbf{G}{=}0,\mathbf{G'}{=}0}(q{\rightarrow}0,\omega)]^{-1}$, is also computed with the BSE.

The tetrahedron method is employed for integration over the Brillouin zone to calculate the optical spectrum. When
calculating the dielectric response within BSE, the valence and conduction states that form the two-particle Hamiltonian
are increased until the two-particle eigenvalues converge within an accuracy of 10 meV. For ReS$_{2}$ the excitons are
essentially Wannier-Mott~\cite{qiu2013optical} in nature and only the states at the valence band top and conduction band
bottom contribute to their formation, so the convergence in the two-particle Hamiltonian size is much faster compared to
the cases of CrX$_{3}$~\cite{acharya2022real} where the excitons have Frenkel character and many valence and conductions
bands over several electron volts form them. However, for the present work, our focus is the plasmonic response, while
we will note the excitonic binding energies in different layered variants later.

Table (\ref{tab:Lattice_Parameters}) contains the lattice parameters for $\text{ReS}_2$ used throughout the
calculations.

\begin{table*}[!ht]
	\centering
	\caption{Lattice parameters of bulk, monolayer (ML) and bilayer (BL) $\text{ReS}_2$}
	\label{tab:Lattice_Parameters}
	\begin{tabular*}{\textwidth}{@{\extracolsep{\stretch{1}}}*{7}{r}@{}}
		\toprule
		\toprule
		Structure & a(\si{\angstrom}) & b(\si{\angstrom}) & c(\si{\angstrom}) & $\alpha$(\textdegree)& $\beta$(\textdegree)& $\gamma$(\textdegree) \\
		\midrule
		Bulk        & 6.41695      & 6.52047    & 7.28252 & 91.8128 & 103.5630 &  118.8390\\

		\\

		ML          & 6.41910      & 6.52306    & 45 & 90.7434 & 95.7909 &  118.8366\\

		\\

		BL          & 6.41910      & 6.52306    & 28.66062  & 84.0127  & 89.7412 &  61.1634\\
		\bottomrule
	\end{tabular*}

\end{table*}

\begin{table}[!ht]
	\centering
	\setlength\tabcolsep{0pt}
	\caption{Bandgap of bulk, BL and ML variants of $\text{ReS}_2$ at different levels of theory (with spin-orbit coupling). The gap increases from LDA to QS\emph{GW} level. The effect of screening is only moderately increased when two particle interactions are added (via a BSE, W $\rightarrow \hat{W}$) are added, thus only weakly decreasing the QS\emph{GW} band gap. }
	\label{tab:results_band_structure}
	\begin{tabular}{ p{2cm} p{2cm} p{2cm} p{2cm}}

		\multicolumn{4}{c}{\textbf{Band Gap(\si{\electronvolt})}}\\
		\toprule

		\textbf{Theory} & \textbf{LDA}  & \textbf{QS\emph{GW}} & \textbf{QSG}$\widehat{\textbf{W}}$ \\
		\toprule
		\toprule
		\textbf{Bulk}          & 1.15     & 1.75 & 1.7  \\

		\\

		\textbf{ML}          & 1.29      & 2.75 & 2.66 \\

		\\

		\textbf{BL}          & 1.23      &  2.35 &  2.3  \\

		\bottomrule
	\end{tabular}
\end{table}

\section{\textsc{Results and Discussions}}
\subsection{\textsc{Quasiparticle Energies and Band Structure}}\label{sec: band-str}

\begin{figure*}[htb]

	\subfloat[\label{bse-mono:a}]{%
		\includegraphics[scale=0.37]{bse-mono.pdf}%
	}
	\subfloat[\label{bse-bi:b}]{%
		\includegraphics[scale=0.37]{bse-bi.pdf}%
	}
	\subfloat[\label{bse-bulk:c}]{%
		\includegraphics[scale=0.37]{bse-bulk.pdf}%
	}

	\caption{\bse band structures (with spin-orbit coupling) with contributions from Re (Red) and S (green). The nature of the band gap is indirect for all the variants of $\text{ReS}_2$ with values,  \SI{2.66}{\electronvolt} for ML (left), \SI{2.3}{\electronvolt} for BL (center) and \SI{1.7}{\electronvolt}  for bulk (right).}
	\label{fig:bse-bandgap}

\end{figure*}

\begin{table*}[!ht]
	\centering
	\setlength\tabcolsep{0pt}
	\caption{Dielectric Constant(real part of the dielectric response at $\omega =0$) calculated using \bse. As the
          dimensionality of $\text{ReS}_2$ is lowered, the frequency range where Re($\epsilon_{xx}$) $\cdot$
          Re($\epsilon_{yy}$) is less than zero becomes smaller.}
	\label{tab:optical-spectra}
	\begin{tabular}{ p{1.5cm} p{1.5cm} p{1.5cm} p{1.5cm} p{3.5cm} p{2.5cm}}

		\multicolumn{4}{c}{\textbf{\thead{Dielectric Constant \\(\si{\electronvolt})}}}\\
		\toprule

	 & $\epsilon_{\infty}^{xx}$  & $\epsilon_{\infty}^{yy}$ & $\epsilon_{\infty}^{zz}$ & \textbf{\thead{Plasmonic frequency \\ range \\(\si{\electronvolt})}} & \textbf{\thead{Exciton \\  Binding energy \\(\si{\electronvolt}) }}\\
		\toprule
		\toprule
		\textbf{Bulk}          & 9.67     & 9.37 & 6.21 & 6.02 - 6.78 (0.76) & -  \\

		\\

		\textbf{BL}          & 6.97      &  7.09 &  2.66 & 6.65 - 7.08 (0.43) & 0.3 \\

		\\

		\textbf{ML}          & 2.97      & 3.19 & 1.42 & - & 0.74 \\

		\bottomrule
	\end{tabular}
\end{table*}

\begin{figure*}[htb]

	\subfloat[\label{mono-xx-yy-bse-real:a}]{%
		\includegraphics[scale=0.33]{mono-xx-yy-bse-real.pdf}%
	}
	\subfloat[\label{bi-xx-yy-bse-real:b}]{%
		\includegraphics[scale=0.33]{bi-xx-yy-bse-real.pdf}%
	}
	\subfloat[\label{bulk-xx-yy-bse-real:c}]{%
		\includegraphics[scale=0.33]{bulk-xx-yy-bse-real.pdf}%
	}

	\subfloat[\label{mono-xx-yy-bse-img:d}]{%
		\includegraphics[scale=0.33]{mono-xx-yy-bse-img.pdf}%
	}
	\subfloat[\label{bi-xx-yy-bse-img:e}]{%
		\includegraphics[scale=0.33]{bi-xx-yy-bse-img.pdf}%
	}
	\subfloat[\label{bulk-xx-yy-bse-img:f}]{%
		\includegraphics[scale=0.33]{bulk-xx-yy-bse-img.pdf}%
	}

	\caption{Real part($\epsilon_1$) of the dielectric response (top row) and imaginary part (bottom row) along $x$ and $y$ direction for ML (left), BL (center) and bulk (right). The vertical dashed line marks the plasmonic frequency range, where $\epsilon_1^{xx}(\omega) \times \epsilon_1^{yy}(\omega) < 0$.}
	\label{fig:bse-real-img}

\end{figure*}

The nature of the bandgap of $\text{ReS}_2$ has been widely debated in the literature. In a typical TMD family, the bandgap is direct when their thickness is reduced to the monolayer, ensuring that coupling with light is strong.  One study from 2014 \cite{tongay2014monolayer} reported direct bandgap for bulk ReS$_{2}$, thus generating considerable interest in the system. However, both older studies such as \cite{ho1997optical,ho1997temperature} and more recent studies such as \cite{aslan2016linearly, gutierrez2016electroluminescence} report that the bulk $\text{ReS}_2$ is an indirect-bandgap semiconductor.

The electronic band structures (with spin-orbit coupling included) for bulk, BL and ML are shown in
Fig. \ref{fig:bse-bandgap} and the band gaps at different levels of theory are summarized in Table
\ref{tab:results_band_structure}. The free-standing ML of $\text{ReS}_2$ has been simulated for these calculations with
the parameters shown in Table \ref{tab:Lattice_Parameters}. Similar to prior work on monolayers of chromium trihalides
\cite{acharya2021electronic}, we check for convergence and scaling of band gap and the dielectric constant
$\epsilon_{\infty}$ with vacuum size. We obtain the LDA band gap of $\sim$\SI{1.29}{\electronvolt} which is
significantly lower than the QS\emph{GW} band gap of $\sim$\SI{2.75}{\electronvolt}. LDA is known to underestimate the band gaps in semiconductors, and the enhancement in
the QS\emph{GW} band gap relative to LDA is standard~\cite{van2006quasiparticle}.  QS\emph{GW} usually
  overcorrects the gap because \emph{W} is universally too large within the random phase approximation (RPA),
and for the same reason it underestimates the dielectric constant $\epsilon_{\infty}$
  \cite{cunningham2018effect}.  Adding ladders largely eliminates both tendencies.  In the present case extending
$\text{QS}GW{\rightarrow}\text{QS}G\widehat{W}$ causes only a modest reduction in the gap, to
$\sim$\SI{2.66}{\electronvolt}, suggesting insignificant corrections to the self-energy originating from BSE. This
result is consistent with previous theoretical work \cite{zhong2015quasiparticle} on ML-$\text{ReS}_2$. The self-energy
and the reduced screening increase the band gap and modify the band topology, which is observed in the
ML-$\text{ReS}_2$. The band gap at the level of LDA is direct, but the two-particle interactions lower the valence band
maxima (VBM), which was at the $\Gamma$ point by about \SI{150}{\milli \electronvolt}. A similar kind of change in band
topology has been observed for $\text{ReSe}_2$ in \cite{zhong2015quasiparticle} but is absent in other some theoretical
work \cite{echeverry2018theoretical}.

We obtain a similar band topology in bulk and BL variants of $\text{ReS}_2$; however, the band gap values are
$\sim$\SI{1.7}{\electronvolt} and $\sim$\SI{2.3}{\electronvolt} respectively. The nature of band gap in $\text{ReS}_2$
is different from more commonly studied semiconducting TMDs (e.g., MoS$_{2}$, WS$_{2}$, etc.), where the bulk and
few-layer variants show indirect band gap and ML is direct. In this work, we observe a direct band gap at the LDA level
and an indirect band gap for all the variants at the \bse level. This nature of the quasiparticle band gap does not
conflict with experimental measurement because, different from our calculated free-standing cases, these measured
samples are on substrates and are inevitably doped. The resulting self-energy corrections will be reduced under these
conditions, resulting in a slightly indirect band gap for the measured samples. We note that the direct to indirect
transition may not be sharp because the energy difference between the direct and indirect gap is small and external
perturbations can affect the conclusion. However, our observations on the indirect nature of the band gap is important for many reasons. In most monolayer TMDs, the gap is direct and it makes them sufficiently bright and also the electron-hole radiative lifetimes are extremely small (often in picoseconds). While in systems with indirect band gaps, radiative and non-radiative processes compete, since the indirect states have longer lifetimes making them candidates for optoelectronics and photovoltaics~\cite{pnas}.

\subsection{\textsc{Optical Absorption Spectra : Hyperbolic Plasmons}}\label{sec: optical-spectra}

\begin{figure*}[htb]

	\subfloat[\label{mono-plasmon:a}]{%
		\includegraphics[scale=0.33]{mono-plasmon.pdf}%
	}
	\subfloat[\label{bi-plasmon:b}]{%
		\includegraphics[scale=0.33]{bi-plasmon.pdf}%
	}
	\subfloat[\label{bulk-plasmon:c}]{%
		\includegraphics[scale=0.33]{bulk-plasmon.pdf}%
	}

	\caption{The product of the real part of the dielectric response along X($\epsilon_1^{xx}(\omega)$) and Y($\epsilon_1^{yy}(\omega)$) direction for ML(left),BL(center) and Bulk(right). The vertical dashed line marks the frequency range in which the product is negative. This frequency range is hyperbolic plasmonic frequency range.}
	\label{fig:plasmons}

\end{figure*}

The anisotropic optical absorption has been previously studied using DFT \cite{ho2004optical}. Experimentally
\cite{aslan2016linearly} it is demonstrated that the reduced crystal symmetry of $\text{ReS}_2$ leads to anisotropic
optical properties that persist from the bulk down to the monolayer limit. The absence of excitonic correlations and
underestimated band gaps in LDA studies hide several physical consequences in $\text{ReS}_2$. Advanced theoretical
studies such as \cite{zhong2015quasiparticle} tackle anisotropic optical responses at the BSE level for monolayer
$\text{ReS}_2$, where ladder vertex corrected optical properties are computed on top of a single shot DFT based
G$_{0}$W$_{0}$ one-particle description.

The large structural anisotropy in 2D materials, for example, a 4:3 anisotropy of the in-plane lattice constants in black
phosphorus~\cite{nemilentsau2016anisotropic} and solid
nitrogen~\cite{rudenko2022electronic,ji2020nitrogen,PhysRevLett.124.216001}) makes them perfect candidates for
hyperbolic materials and a natural place to look for HP. This offers new possibilities as hyperbolic materials showcase a wide variety of interesting properties, such as modes which transport heat by photon tunnelling with a high efficiency close to the theoretical limit\cite{biehs2012hyperbolic}, and broadband absorption \cite{riley2017near}. We first define the condition for hyperbolicity. The hyperbolic region appears when
\begin{equation}
		\epsilon_1^{xx}(\omega) \cdot \epsilon_1^{yy}(\omega) < 0
\end{equation}
where, $\epsilon_1^{xx}(\omega)$ and $\epsilon_1^{yy}(\omega)$ are the real part of the dielectric response along the x and y direction respectively. We assume that $\epsilon_1^{xy}(\omega)=0$ by symmetry and, thus, $x$ and $y$ are the principal directions of the dielectric permeattivity. For different variants of $\text{ReS}_2$, the real and imaginary part of dielectric response is plotted in Fig. \ref{fig:bse-real-img}. We observe that significant difference in optical response for incident polarized light along different directions. While ML-$\text{ReS}_2$ (Fig. \ref{fig:bse-real-img}, left panel), hosts some strongly bound anisotropic excitons deep inside the one-particle gap, HP are absent. For the BL-$\text{ReS}_2$ (Fig. \ref{fig:bse-real-img}, center panel), the Re($\epsilon_{yy}$) goes to negative at $\omega= \SI{6.65}{\electronvolt}$, which results in a hyperbolic region starting at that frequency with an energy window of \SI{0.43}{\electronvolt}. For bulk-$\text{ReS}_2$ this energy window increases to \SI{0.76}{\electronvolt}. The sign change is key to the appearance of the hyperbolic region and it becomes more apparent in Fig. \ref{fig:plasmons} where we plot the product $\epsilon_1^{xx}(\omega) \times \epsilon_1^{yy}(\omega)$ which becomes negative in the energy window.  $\epsilon_2$ remains large for both bulk and BL in the hyperbolic energy window suggesting large damping of the HP. In the BL, these plasmons are less damped compared to bulk. Also, note that these HP in ReS$_{2}$ are in ultraviolet range, in contrast to the infrared HP in CuS nanocrystals~\cite{cordova2019anisotropic}.

The inherent anisotropy in ReS$_{2}$ provides an opportunity to tune its magnitude by applying strain. We apply uni-directional strains ($\gamma$) along $x$ and $y$ respectively and explore the hyperbolic region. We apply up to 4\% strain and see that the ML never hosts HP. However, in bulk  on application of compressive unidirectional strain ($\gamma_{x}$) along $x$, the HP window increases upto $\sim$1.3 eV. This enhanced HPs also have lesser damping compared to the un-strained compound. On the other hand, $\gamma_{y}$ reduces the HP window. In the BL, application of $\gamma_{x}$, almost entirely kills the HP window while $\gamma_{y}$ enhances the HP window and also hosts less damped plasmonic modes. In short, we observe that while strain can be used to tune the hyperbolic energy window and the lifetimes of the plasmons, we could not produce HPs in the monolayer sample under any condition. However, the HPs remain pretty robust in both the bulk and BL variants, also their stability could be enhanced on selective applications of uni-directional strain.

\section{Conclusion} \label{sec: conclusion}

Anisotropy is a key to tuning material properties. Discontinuities at surfaces, residual strains and metamaterials
have been used as platforms for realizing anisotropic optical properties. Naturally occurring structurally anisotropic materials are not necessarily hyperbolic always.

In this work we show that the structural anisotropy in $\text{ReS}_2$, even though much weaker compared to materials like black Phosphorus and solid Nitrogen, leads to the occurrence of hyperbolic plasmons in a narrow energy window. The plasmonic resonances can be tuned by controlling the number of layers in the far ultraviolet frequency range. The ability of such a tunability of the plasmons opens up new opportunities in regard to optoelectronic devices. We further show that the hyperbolic region and its stability can be enhanced by unidirectional strains.

\section{Acknowledgement} \label{sec: Acknowledgement}
MIK and SA are supported by the ERC Synergy Grant, project 854843 FASTCORR (Ultrafast dynamics of correlated
electrons in solids). MvS and DP were supported by the Computational Chemical Sciences program
within the Office of Basic Energy Sciences, U.S. DOE under Contract
No. DE-AC36-08GO28308. This research used resources of the National
Energy Research Scientific Computing Center (NERSC), award
BES-ERCAP0021783, under DOE Contract No. DE-AC02-05CH11231. We acknowledge PRACE
for awarding us access to Irene-Rome hosted by TGCC, France and Juwels Booster and Cluster, Germany.