††thanks: These authors contribute to this work equally.††thanks: These
authors contribute to this work equally.

# Precise Fermi level engineering in a topological Weyl semimetal via fast ion
implantation

Manasi Mandal  Quantum Measurement Group, MIT, Cambridge, MA 02139, USA
Department of Nuclear Science and Engineering, MIT, Cambridge, MA 02139, USA
Abhijatmedhi Chotrattanapituk  Quantum Measurement Group, MIT, Cambridge, MA
02139, USA  Department of Electrical Engineering and Computer Science, MIT,
Cambridge, MA 02139, USA  Kevin Woller  Department of Nuclear Science and
Engineering, MIT, Cambridge, MA 02139, USA  Haowei Xu  Department of Nuclear
Science and Engineering, MIT, Cambridge, MA 02139, USA  Nannan Mao  Department
of Electrical Engineering and Computer Science, MIT, Cambridge, MA 02139, USA
Ryotaro Okabe  Quantum Measurement Group, MIT, Cambridge, MA 02139, USA
Department of Chemistry, MIT, Cambridge, MA 02139, USA  Artittaya Boonkird
Quantum Measurement Group, MIT, Cambridge, MA 02139, USA  Department of
Nuclear Science and Engineering, MIT, Cambridge, MA 02139, USA  Thanh Nguyen
Quantum Measurement Group, MIT, Cambridge, MA 02139, USA  Department of
Nuclear Science and Engineering, MIT, Cambridge, MA 02139, USA  Nathan C.
Drucker  Quantum Measurement Group, MIT, Cambridge, MA 02139, USA  School of
Engineering and Applied Sciences, Harvard University, Cambridge, MA 02138, USA
Takashi Momiki  Aisin Technical Center of America, Inc, CA 95110, USA  Ju Li
Department of Nuclear Science and Engineering, MIT, Cambridge, MA 02139, USA
Department of Materials Science and Engineering, MIT, Cambridge, MA 02139, USA
Jing Kong  Department of Electrical Engineering and Computer Science, MIT,
Cambridge, MA 02139, USA  Mingda Li††<EMAIL_ADDRESS><EMAIL_ADDRESS>Quantum
Measurement Group, MIT, Cambridge, MA 02139, USA  Department of Nuclear
Science and Engineering, MIT, Cambridge, MA 02139, USA

###### Abstract

The precise controllability of the Fermi level is a critical aspect of quantum
materials. For topological Weyl semimetals, there is a pressing need to fine-
tune the Fermi level to the Weyl nodes and unlock exotic electronic and
optoelectronic effects associated with the divergent Berry curvature. However,
in contrast to 2D materials, where the Fermi level can be controlled through
various techniques, the situation for bulk crystals beyond laborious chemical
doping poses significant challenges. Here, we report the meV-level ultra-fine-
tuning of the Fermi level of bulk topological Weyl semimetal TaP using
accelerator-based high-energy hydrogen implantation and theory-driven
planning. By calculating the desired carrier density and controlling the
accelerator profiles, the Fermi level can be fine-tuned from 5 meV to only
$\sim$0.5 meV (DFT calculations) away from the Weyl nodes. The Weyl nodes are
preserved, while the carrier mobility is largely retained. Our work
demonstrates the viability of this generic approach to tune the Fermi level in
semimetal systems and could serve to achieve property fine-tuning for other
bulk quantum materials with ultrahigh precision.

The precise controllability of the Fermi level and carrier doping in a solid-
state material holds paramount significance in condensed matter physics and
materials science. In microelectronics such as high-electron-mobility
transistors, the ability to fine-tune the Fermi level is crucial to control
the operation thresholds and high-frequency responses [1, 2, 3]. In strongly
correlated systems such as high-temperature cuprate superconductors, achieving
optimal hole doping is essential for enhancing the superconducting critical
temperature [4, 5, 6]. As for energy harvesting applications such as
thermoelectrics, a well-adjusted Fermi level strikes a delicate balance
between electrical conductivity and thermopower, which leads to an optimized
power factor [7, 8, 9, 10]. A particularly noteworthy utilization for Fermi
level fine-tuning lies in topological materials, where only a precisely
positioned Fermi level can unveil the topological properties. For instance, in
topological insulators, dissipationless electronic states such as quantum
anomalous Hall states can only manifest when the Fermi level is located
precisely within the small surface bandgap. This can only be achieved through
a strategic blend of chemical doping [11, 12, 13] and electrostatic gating in
thin heterostructures [14, 15, 16, 17]. For topological Weyl semimetals (WSM),
exotic phenomena like large charge-to-spin interconversion [18, 19], strong
higher-order photoresponse [20, 21], anomalous Nernst effect [22, 23], and
quantized thermoelectric Hall effect [24, 25] emerge due to a divergent Berry
curvature, which can only happen when the Fermi level is tuned to be exactly
at the Weyl nodes.

However, unlike in lower dimensional materials, in which the Fermi level can
be continuously fine-tuned through a variety of accessible techniques like
electrostatic and ionic gating, achieving this level of control in three-
dimensional bulk materials has presented a longstanding challenge. Chemical
doping can only be performed during the synthesis and offers limited
precision, which does not meet the stringent requirement for topological
semimetals. The electronic states in topological semimetals, such as from the
Weyl and Dirac bands, have a linear dispersion with zero density-of-states at
the topological singularities, making the Fermi level location extremely
sensitive to the carrier concentration [26, 27]. Therefore, there is an urgent
need to attain enhanced fine-tunability of topological semimetals for more
practical energy and information applications.

abcd Figure 1: a) Tandem accelerator schematic. It includes an ion source to
produce negative ions. Selected ion species are injected toward the terminal
using a low-energy magnet. The ion beam is focused by a magnetic quadrupole
lens and the desired ion species are directed into the beamline for the
experiment. b) Crystal structure of the Weyl semimetal TaP. c) Fast Fourier
transform of the longitudinal magnetoresistance of TaP reveals quantum
oscillations with specific frequencies (reproduced from [24]). d) Illustration
of the irradiation effect on the Fermi level of TaP. The red arrow indicates
the shift of the Fermi level towards a Weyl node.

In this work, we use high-energy, accelerator-based hydrogen implantation to
achieve ultra-fine carrier doping in WSM of tantalum phosphide (TaP) [24].
Although hydrogenation and low-energy ion implantation have been used to tune
materials properties [28, 29, 30, 31, 32, 33], there are a few key distinct
features presented in the current setup. Foremost, the fine controllability of
the advanced accelerator technique enables an ultra-fine tuning of the Fermi
level down to the meV regime, far beyond the previous approaches. In addition,
the high energy beam also enables the doping of a bulk crystal beyond a thin
film with tens of nanometers [30, 31, 34]. Moreover, thanks to the rapid
development of the accelerator technology with high precision control of the
accelerator beam energy and flux of ions, we were able to achieve a priori
doping planning through a combination of DFT and Monte Carlo calculations,
which eventually shows quantitative agreement with actual experiments. This
sets up an example for ultrafine doping planning experiments. The Weyl nodes
are preserved after doping, as shown from the quantum oscillation
measurements, while the temperature to reach the charge neutral point is
nearly doubled, without lowering much on the carrier mobility.

abTa2P2H in Ta16P16 Figure 2: Ab initio calculation of Weyl nodes’ energy
relative to Fermi level ($E_{f}$) of a) pristine TaP, and b) H doped TaP
($2\times 2\times 2$ TaP supercell with one H insertion). Each plot point
represents a Weyl node placed in its energy ascending order from left to
right. The inset in b shows the zoom-in view of the Weyl nodes with energy
near Fermi level

The doping experiment uses the Tandem accelerator located at the Center for
Science and Technology with Accelerators and Radiation (CSTAR). The schematic
used in our experiment is shown in Fig. 1a. This type of accelerator consists
of two successive linear accelerators with the same power supply. A cesium
sputtering source is used as the ion source to generate negative ions from
most of the elements with low vapor pressure (excluding noble gases). This
source works by cesium sputtering a solid target, which produces negative ions
that are then subjected to electron exchange with neutral cesium. Negative
ions are initially introduced into the first stage (blue arrow), where they
are accelerated toward the positive high-voltage (HV) terminal. Subsequently,
depending on the experimental preferences, nitrogen gas can be injected into
the HV terminal to strip some of the electrons from the ions, which transforms
them into positive ions. As a result, the Tandem accelerator enables a
remarkable capability to dope tens of different elements in the period table
of elements. We have chosen negative hydrogen as the dopant, which offers
several advantages in materials science due to its small atomic size, i.e., it
does not significantly alter the structure of pristine materials [35].
Electrochemical hydrogenation can also be used for electron doping of
semiconductors with good doping regime controllability, but not high-precision
doping level control [36]. In addition, hydrogen can exhibit a range of exotic
properties as a dopant that makes it particularly interesting for research and
applications, including metal-to-insulator transitions [37], magnetostriction
[38], enhanced superconductivity [39], and paraelectric-to-ferroelectric phase
transitions [40].

abcdeS3mS3mS3mS3mS3m

fghijS20mS20mS20mS20mS20m

Figure 3: a) Variation of magnetoresistance (MR) with magnetic field ($B$) at
different temperatures for sample S3m. b) Background-subtracted MR
($\Delta$MR) of a versus 1/B for a few temperatures that exhibit Shubnikov-de
Haas (SdH) oscillations. This $\Delta$MR is calculated using q polynomial-
fitted background subtraction. c) Fast Fourier transform (FFT) of b. d)
Example of a Landau level fan diagram analysis which indicates the presence of
lower frequencies in b at 5 K than TaP [24]. e) Scatter plot of the
oscillation frequencies extracted from b using various background subtraction
methods at 5 K. f)-j) are the same as a)-e), but for the S20m sample.

The earlier study [24] highlights TaP as a remarkable topological Weyl
semimetal with an inversion symmetry-breaking crystal structure (Fig. 1b). It
exhibits a non-saturating thermopower, a giant magnetoresistance, and
quantized thermoelectric Hall effect. Shubnikov-de Haas oscillations (SdH) are
observed in the background-subtracted MR ($\Delta$MR) data at temperatures
below 25 K. The Fourier transform of $\Delta$MR showed four pockets, with two
small carrier pockets, one with a low-frequency Fα = 4 T and the other with Fβ
= 18 T indicating the presence of the W2 Weyl node, in addition to the
electron pocket from the W1 Weyl node contributing to Fβ (Fig. 1c) [24].
Introducing H- ions through irradiation is expected to shift the Fermi level
closer to the Weyl node by carrier doping, as illustrated in the schematic
diagram in Fig. 1d.

By employing the Vienna ab initio simulation package (VASP), we conducted
density functional theory (DFT) calculations. These calculations provide a
qualitative insight into the energy redistribution of Weyl node energy levels
upon H doping. The energy of the Weyl nodes with respect to the Fermi level
($E_{f}$) is shown in Fig. 2a for TaP and Fig. 2b for the doped TaP system
(one H in Ta16P16). Here, the horizontal axis is the index of the Weyl nodes.
The Weyl nodes are ranked according to their energy, so the index itself does
not carry additional information. Notably, the doped TaP system exhibits Weyl
nodes that are distributed somewhat continuously. Comparing the doped TaP to
the pristine TaP, our calculations indicate that certain Weyl nodes in the
doped system lie merely 0.5 meV away from the Fermi level, whereas, in the
pristine system, they are much further at 10 meV from the Fermi level. A
theoretical model and ab initio calculations estimating the correlation
between hydrogen doping concentration in TaP and the shift in the Fermi level
are discussed in Supplementary Information (SI).

To calculate the implanted ion energy, penetration depth within the crystal,
and damage to the target material, the stopping and range of ions in matter
software (SRIM/TRIM) [41] was employed. The detailed calculations are
presented in the SI. Single crystals of TaP were irradiated with 20 keV H-
ions with varying durations to achieve the desired doping concentration. The
summary of different samples with their respective dose of irradiation and
time duration is presented in Table 1.

The irradiated samples named S3m and S20m exhibit a significantly high
magnetoresistance (MR) on the order of 104% and accompanied by SdH
oscillations at temperatures below 15 K (see Fig. 3a, and Fig. 3f,
respectively). In contrast, the S2h sample does not display such oscillations
even up to 9 T (details in SI). It is important to note that the analysis of
carrier pockets based on quantum oscillations can be influenced by the choice
of the background for MR. Therefore, we calculated the background-subtracted
MR ($\Delta$MR) using different background subtraction methods (details in
SI).

Sample | Dose of irradiation | Irradiation time
---|---|---
| (cm-2) | (minutes)
S3m | 1.25$\times$1015 | 3
S20m | 1.23$\times$1016 | 20
S2h | 1.23$\times$1017 | 120
Table 1: Summary of TaP samples irradiated with H- ions for different
durations. abPristineS3mS20mPristineS3mS20m Figure 4: a) Temperature-dependent
profiles of electron ($n_{e}$) and hole concentrations ($n_{h}$). b)
Temperature-dependent trends of electron ($\mu_{e}$) and hole mobilities
($\mu_{h}$). The top profiles serve as references for TaP. Our analysis
involves the extraction of equivalent profiles from the two-band model
fitting, applied to the $S_{3m}$ and $S_{20m}$ samples as annotated in the
plots.

Fig. 3b, and Fig. 3g illustrate the $\Delta$MR obtained with a polynomial-
fitted background subtraction method, while Fig. 3c, and Fig. 3h show the
corresponding fast Fourier transform (FFT) results for S3m and S20m,
respectively. These FFT plots reveal a noticeably lower oscillation frequency
compared to the pristine TaP [24]. Furthermore, Fig. 3d, and Fig. 3i present
the results of the Landau level (LL) fan diagram analysis which indicates the
presence of more oscillations with lower frequencies. The oscillation
frequencies have been repeatedly extracted using various background
subtraction methods, and scatter plotted together to show agreement between
methods (see Fig. 3e, and Fig. 3j). This finding supports and validates our
DFT calculations concerning the continuous energy distribution of Weyl nodes
in proximity to the Fermi surface compared to pristine TaP.

To find the effect of irradiation on the carrier concentration and the
mobility Fig. 4, the longitudinal and transverse conductivity data were fitted
with a two-band model. The details of the data extraction from these
measurements are discussed in the SI. The compensation temperature increased
to 160 K and 155 K for S3m and S20m, respectively, whereas it is 110 K for
pristine TaP [24]. Fig. 5 shows the normalized carrier concentration $n_{T}$
extracted from the simple extended Kohler’s rule. The consistent temperature-
dependent trend observed in our analysis, in accordance with the fitted two-
band model, serves as compelling evidence for the accuracy and validity of our
findings (see SI for detailed calculation).

In summary, we have addressed a novel method for ultra-fine Fermi level tuning
in the prototype bulk WSM TaP using high-energy hydrogen implantation
facilitated by accelerator-based techniques. Our transport measurements
demonstrate a successful increment in the charge neutral point temperature
through precise tuning of the Fermi level in proximity to the Weyl nodes.
These experimental results are achievable by the guidance from DFT
calculations. The approach presented in this study represents a highly
controllable and universally applicable technique for fine-tuning the Fermi
level and quantum orderings in bulk quantum materials. By enabling precise
control over the Fermi level in bulk crystals, the work could open up new
possibilities for exploring and manipulating the unique properties of
topological WSM and other quantum materials. This advancement holds
significant promise for the development of future quantum technologies and
applications.

abS3mS20m Figure 5: Direct comparison of the temperature dependence of the
normalized adjusted parameter $n_{T}$, as extracted from both the extended
Kohler’s rule and the two-band model for a) S3m, and b) S20m. The plots show
good agreement between simple extended Kohler’s rule and fitted two-band
model.

###### Acknowledgements.

Acknowledgements
The authors Manasi Mandal and Abhijatmedhi Chotrattanapituk contributed
equally to this work. MM and AC acknowledge support from the US Department of
Energy (DOE), Office of Science (SC), Basic Energy Sciences (BES), Award No.
DE-SC0020148. TN is supported by NSF Designing Materials to Revolutionize and
Engineer our Future (DMREF) Program with Award No. DMR-2118448. TB is
supported by the NSF Convergence Accelerator Award No. 2235945. JL
acknowledges support by NSF, DMR-2132647. ML acknowledges the support from
Aisin Inc, the Class of 1947 Career Development Professor Chair, and the
support from R. Wachnik.

## References

  * Hemaja and Panda [2022] V. Hemaja and D. Panda, Silicon 14, 1873 (2022).
  * Lee  _et al._ [2016] H. Lee, J. Perozek, L. D. Rosario, and C. Bayram, Sci. Rep. 6, 37588 (2016).
  * Khan  _et al._ [1993] M. A. Khan, A. Bhattarai, J. N. Kuznia, and D. T. Olson, Appl. Phys. Lett. 63, 1214 (1993).
  * Hossain  _et al._ [2008] M. A. Hossain, J. D. Mottershead, D. Fournier, A. Bostwick, J. L. McChesney, E. Rotenberg, R. Liang, W. N. Hardy, G. A. Sawatzky, I. S. Elfimov, D. A. Bonn, and A. Damascelli, Nat. Phys. 4, 527 (2008).
  * Michon  _et al._ [2023] B. Michon, C. Berthod, C. W. Rischau, A. Ataei, L. Chen, S. Komiya, S. Ono, L. Taillefer, and A. Georges, Nat. Commun. 14, 3033 (2023).
  * Weber  _et al._ [2010] C. Weber, K. Haule, and G. Kotliar, Nat. Phys. 6, 574 (2010).
  * Kim  _et al._ [2013] G. Kim, L. Shao, K. Zhang, and K. P. Pipe, Nat. Mater. 12, 719 (2013).
  * Wu  _et al._ [2017] L. Wu, X. Li, S. Wang, T. Zhang, J. Yang, W. Zhang, L. Chen, and J. Yang, NPG Asia Mater. 9, e343 (2017).
  * Tsuda  _et al._ [2023] S. Tsuda, A. Yoshinari, S. Takezawa, K. Ohishi, N. Nagamura, W. Zhang, Y. Iwasaki, and Y. Takagiwa, Mater. Res. Express 10, 055506 (2023).
  * Pei  _et al._ [2013] Y. Pei, J. He, J. Li, F. Li, Q. Liu, W. Pan, C. Barreteau, D. Berardan, N. Dragoe, and L. Zhao, NPG Asia Mater. 5, e47 (2013).
  * Chang  _et al._ [2013] C.-Z. Chang, J. Zhang, X. Feng, J. Shen, Z. Zhang, M. Guo, K. Li, Y. Ou, P. Wei, L.-L. Wang, Z.-Q. Ji, Y. Feng, S. Ji, X. Chen, J. Jia, X. Dai, Z. Fang, S.-C. Zhang, K. He, Y. Wang, L. Lu, X.-C. Ma, and Q.-K. Xue, Science 340, 167 (2013).
  * Kang  _et al._ [2023] S. Kang, S. Kang, H. Kim, and J. Yu, npj 2D Mater. Appl. 7, 13 (2023).
  * Zhang  _et al._ [2013] X. Zhang, L. Liu, and W. Liu, Sci. Rep. 3, 2908 (2013).
  * Ferguson  _et al._ [2023] G. M. Ferguson, R. Xiao, A. R. Richardella, D. Low, N. Samarth, and K. C. Nowack, Nat. Mater. 22, 1100 (2023).
  * Zhang  _et al._ [2023] C. Zhang, T. Zhu, T. Soejima, S. Kahn, K. Watanabe, T. Taniguchi, A. Zettl, F. Wang, M. P. Zaletel, and M. F. Crommie, Nat. Commun. 14, 3595 (2023).
  * Deng  _et al._ [2020] Y. Deng, Y. Yu, M. Z. Shi, Z. Guo, Z. Xu, J. Wang, X. H. Chen, and Y. Zhang, Science 367, 895 (2020).
  * Lin  _et al._ [2022] W. Lin, Y. Feng, Y. Wang, J. Zhu, Z. Lian, H. Zhang, H. Li, Y. Wu, C. Liu, Y. Wang, J. Zhang, Y. Wang, C.-Z. Chen, X. Zhou, and J. Shen, Nat. Commun. 13, 7714 (2022).
  * Zhao  _et al._ [2020] B. Zhao, D. Khokhriakov, Y. Zhang, H. Fu, B. Karpiak, A. M. Hoque, X. Xu, Y. Jiang, B. Yan, and S. P. Dash, Phys. Rev. Research 2, 013286 (2020).
  * Mendes  _et al._ [2022] J. B. S. Mendes, A. S. Vieira, R. O. Cunha, S. O. Ferreira, R. D. Dos Reis, M. Schmidt, M. Nicklas, S. M. Rezende, and A. Azevedo, Adv. Mater. Interfaces 9, 2201716 (2022).
  * Osterhoudt  _et al._ [2019] G. B. Osterhoudt, L. K. Diebel, M. J. Gray, X. Yang, J. Stanco, X. Huang, B. Shen, N. Ni, P. J. W. Moll, Y. Ran, and K. S. Burch, Nat. Mater. 18, 471 (2019).
  * Ma  _et al._ [2019] J. Ma, Q. Gu, Y. Liu, J. Lai, P. Yu, X. Zhuo, Z. Liu, J.-H. Chen, J. Feng, and D. Sun, Nat. Mater. 18, 476 (2019).
  * Li  _et al._ [2020] P. Li, J. Koo, W. Ning, J. Li, L. Miao, L. Min, Y. Zhu, Y. Wang, N. Alem, C. Liu, Z. Mao, and B. Yan, Nat. Commun. 11, 3476 (2020).
  * Roychowdhury  _et al._ [2023] S. Roychowdhury, M. Yao, K. Samanta, S. Bae, D. Chen, S. Ju, A. Raghavan, N. Kumar, P. Constantinou, S. N. Guin, N. C. Plumb, M. Romanelli, H. Borrmann, M. G. Vergniory, V. N. Strocov, V. Madhavan, C. Shekhar, and C. Felser, Adv. Sci. 10, 2207121 (2023).
  * Han  _et al._ [2020] F. Han, N. Andrejevic, T. Nguyen, V. Kozii, Q. Nguyen, T. Hogan, D. Z., P. S., S. B., A. A., A. E., C. S., H. S., F. L., and L. M., Nat. Commun. 11, 6167 (2020).
  * Scott  _et al._ [2023] E. F. Scott, K. A. Schlaak, P. Chakraborty, C. Fu, S. N. Guin, S. Khodabakhsh, A. E. Paz y Puente, C. Felser, B. Skinner, and S. J. Watzman, Phys. Rev. B 107, 115108 (2023).
  * lohani  _et al._ [2023] H. lohani, P. Foulquier, P. Le Fèvre, F. m. c. Bertran, D. Colson, A. Forget, and V. Brouet, Phys. Rev. B 107, 245119 (2023).
  * Liu  _et al._ [2023] J. Liu, L. X. Linchao Ding, X. Li, K. Behnia, and Z. Zhu, J. Phys. Condens. Matter 35, 375501 (2023).
  * Xie  _et al._ [2020] L. Xie, Q. Zhu, G. Zhang, K. Ye, C. Zou, O. V. Prezhdo, Z. Wang, Y. Luo, and J. Jiang, J. Am. Chem. Soc. 142, 4136 (2020).
  * Kang  _et al._ [2022] S. Kang, S. Jang, A. N. Morozovska, O. Kwon, Y. Jin, H. Kim, H. Bae, C. Wang, H. Yang, A. Belianinov, S. Randolph, E. A. Eliseev, L. Collins, Y. Park, S. Jo, H. Jung, J. Go, H. W. Cho, Y. Choi, and Y. Kim, Science 376, 731 (2022).
  * Wylie  _et al._ [2022] A. P. C. Wylie, K. B. Woller, S. A. A. A. Dajani, B. R. Dacus, E. J. Pickering, M. Preuss, and M. P. Short, J. Appl. Phys. 132, 045102 (2022).
  * Ferry  _et al._ [2019] S. E. Ferry, C. A. Dennett, K. B. Woller, and M. P. Short, J. Nucl. Mater. 523, 378 (2019).
  * Wang  _et al._ [2023] J. Wang, K. B. Woller, and B. Yildiz, ECS Trans. 111, 809 (2023).
  * Chalifoux  _et al._ [2019] B. D. Chalifoux, Y. Yao, K. B. Woller, R. K. Heilmann, and M. L. Schattenburg, Opt. Express 27, 11182 (2019).
  * Woller [2017] K. B. Woller,  _Characterization of the dynamic formation of nano-tendril surface morphology on tungsten while exposed to helium plasma_ , Ph.D. thesis, Massachusetts Institute of Technology, Department of Nuclear Science and Engineering (2017).
  * Chen  _et al._ [2017] S. Chen, Z. Wang, L. Fan, Y. Chen, H. Ren, H. Ji, D. Natelson, Y. Huang, J. Jiang, and C. Zou, Phys. Rev. B 96, 125130 (2017).
  * Yajima  _et al._ [2018] T. Yajima, G. Oike, S. Yamaguchi, S. Miyoshi, T. Nishimura, and A. Toriumi, AIP Advances 8, 115133 (2018).
  * Wang  _et al._ [2020] Z. Wang, X. Wang, E. Sharman, X. Li, L. Yang, G. Zhang, and J. Jiang, J. Phys. Chem. Lett. 11, 1075 (2020).
  * Nakayama  _et al._ [2018] R. Nakayama, M. Maesato, T. Yamamoto, H. Kageyama, T. Terashima, and H. Kitagawa, Chem. Commun. 54, 12439 (2018).
  * Jiao  _et al._ [2022] N. Jiao, H.-D. Liu, L. Yang, Y.-P. Li, M. Zheng, H.-Y. Lu, and P. Zhang, EPL 138, 46002 (2022).
  * Qing-Xun  _et al._ [2009] Z. Qing-Xun, G. Bo, W. Shu-Biao, B. Fang, G. Li, and L. Bao-Ting, Acta Physico - Chimica Sinica 25, 183 (2009).
  * Ziegler [1984] J. F. Ziegler, SRIM.EXE (1984).
  * Hohenberg and Kohn [1964] P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964).
  * Kohn and Sham [1965] W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965).
  * Kresse and Furthmüller [1996] G. Kresse and J. Furthmüller, Comput. Mater. Sci. 6, 15 (1996).
  * Kresse and Furthmüller [1996] G. Kresse and J. Furthmüller, Phys. Rev. B 54, 11169 (1996).
  * Perdew  _et al._ [1996] J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996).
  * Blöchl [1994] P. E. Blöchl, Phys. Rev. B 50, 17953 (1994).
  * Mostofi  _et al._ [2014] A. A. Mostofi, J. R. Yates, G. Pizzi, Y. S. Lee, I. Souza, D. Vanderbilt, and N. Marzari, Comput. Phys. Commun. 185, 2309 (2014).
  * Wu  _et al._ [2018] Q. S. Wu, S. N. Zhang, H. F. Song, M. Troyer, and A. A. Soluyanov, Comput. Phys. Commun. 224, 405 (2018).
  * Xu  _et al._ [2016] N. Xu, H. M. Weng, B. Q. Lv, C. E. Matt, J. Park, F. Bisti, V. N. Strocov, D. Gawryluk, E. Pomjakushina, K. Conder, N. C. Plumb, M. Radovic, G. Autès, O. V. Yazyev, Z. Fang, X. Dai, T. Qian, J. Mesot, H. Ding, and M. Shi, Nat. Commun. 7, 11598 (2016).
  * ab1 [2023] Materials Project (Accessed in 2023).
  * ab2 [2023] Quantum Espressso pseudopotentials (Accessed in 2023).
  * Xu  _et al._ [2021] J. Xu, F. Han, T.-T. Wang, L. R. Thoutam, S. E. Pate, M. Li, X. Zhang, Y.-L. Wang, R. Fotovat, U. Welp, X. Zhou, W.-K. Kwok, D. Y. Chung, M. G. Kanatzidis, and Z. L. Xiao, Phys. Rev. X 11, 041029 (2021).

##  Appendix A Supplementary Information

###  A.1 DFT

The density functional theory (DFT) [42, 43] calculations are performed using
the Vienna ab initio simulation package (VASP) [44, 45]. Generalized gradient
approximation (GGA) in the form of Perdew-Burke-Ernzerhof (PBE) [46] is used
to treat the exchange-correlation interactions. Core and valence electrons are
treated by the projector augmented wave (PAW) method [47] and a plane wave
basis, respectively. A tight-binding Hamiltonian is built from the DFT results
using the Wannier90 package [48]. The tight-binding Hamiltonian is used to
find Weyl nodes with the help of the WannierTools package [49]. In doped TaP,
the concentration of H is around 2$\times$1021 cm-3.

###  A.2 Theoretical Model of H-doped TaP

a

b

c

Figure S1: Angle-resolved photoemission spectroscopy (ARPES) intensity plot
and its curvature plot with calculated bands for TaP passing through a W1, and
b W2 Weyl nodes (recreated from [50]). c The plot of the equation (7) for the
range $E\in[-50,50]$ meV.

The analysis was initiated with some key aspects of TaP to set up appropriate
assumptions for the model. In TaP conventional unit cell, there are, in total,
24 Weyl nodes, of which 16 are W1 type with energy 24 meV above Fermi level
while the remaining 8 are W2 type with energy 40 meV below Fermi level. The
dispersion is approximately linear near the nodes, especially the regions
between the two node types’ energies. Furthermore, since the dispersion
generally does not change abruptly over a small region in momentum space, we
would approximate all Weyl nodes’ dispersions with tilted 4-dimensional cones.
This means the constant energy surface near nodes’ levels is approximately
ellipsoids. To further simplify the model, we considered the experimental
results of [50]. We adapted Fig. S1a for W1 and Fig. S1b for W2 from figures 2
and 3 of the article [50]. The cones are slightly tilted for both W1 and W2.
Furthermore, the dispersion slopes of principal directions are different as
well. However, because of the large uncertainty of slopes one could extract
from the angle-resolved photoemission spectroscopy (ARPES) experimental
results, and the order of doping magnitude is what we are concerned about (We
ignored the red lines in these figures since they are from the ab initio
calculation). The nodes were approximated with upright but distorted cones,
and their dispersion slopes were extracted from the available ARPES figure of
each node type.

With this approximated model, the constant energy (E) surfaces are ellipsoids
with semi-principal axes

$r_{x,y,z}=\left|\boldsymbol{k}-\boldsymbol{k}_{0}\right|=\frac{\left|\mathbf{E}-\mathbf{E}_{0}\right|}{d_{x,y,z}}.$
(1)

Here, $E_{0}$ is the node’s energy level, $\boldsymbol{k}_{0}$ is the node’s
location, and $d$ is the dispersion slope. This implies that, for each node,
the energy increment from the Fermi level ($E_{f}$) to $E$ would lead to an
expansion in the volume of momentum space located below the Fermi level by

$\displaystyle\left(\frac{4}{3}\right)\pi\dfrac{\left(E-E_{0}\right)^{3}}{d_{x}d_{y}d_{z}}-\left(\frac{4}{3}\right)\pi\dfrac{\left(E_{f}-E_{0}\right)^{3}}{d_{x}d_{y}d_{z}}=\frac{4\pi}{3d_{x}d_{y}d_{z}}\left(\left(E-E_{0}\right)^{3}-\left(E_{f}-E_{0}\right)^{3}\right).$
(2)

The absolute value is omitted from energy differences with $E_{0}$, resulting
in negative volumes for energies below $E_{0}$. This adjustment ensures the
accurate representation of momentum space volume increments. Given the
standard practice of expressing energy relative to the Fermi level, the
formulas were subsequently simplified using this notation as follows:

$\displaystyle\left(\frac{4\pi}{3d_{x}d_{y}d_{z}}\right)\left((E-E_{0})^{3}+E_{0}^{3}\right).$
(3)

Consider energy levels $E_{1}$ and $E_{2}$ corresponding to the nodes of W1
and W2, respectively, where $E_{1}=24$ meV and $E_{2}=-40$ meV. Assuming that
$E$ is close to 0 ($E_{f}$) enough such that only the nodes contribute to the
total momentum space volume increment. The total increment in momentum space
volume (sixteen from W1 and eight from W2) would be

$\displaystyle\left(\frac{32\pi}{3}\right)\left(\frac{2}{{d_{1,x}d_{1,y}d_{1,z}}}\left((E-E_{1})^{3}+E_{1}^{3}\right)+\frac{1}{{d_{2,x}d_{2,y}d_{2,z}}}\left((E-E_{2})^{3}+E_{2}^{3}\right)\right),$
(4)

where $d_{1}$ and $d_{2}$ are average dispersion slopes of W1 and W2 nodes
respectively. This signifies that the number of electronic states increases by

$\displaystyle
2\left(\frac{V}{{8\pi^{3}}}\cdot\frac{32\pi}{3}\right)\left(\frac{2}{{d_{1,x}d_{1,y}d_{1,z}}}\left((E-E_{1})^{3}+E_{1}^{3}\right)+\frac{1}{{d_{2,x}d_{2,y}d_{2,z}}}\left((E-E_{2})^{3}+E_{2}^{3}\right)\right).$
(5)

In an alternative perspective, the requirement arises for TaP to undergo
doping with an additional electron concentration of

$\displaystyle\left(\frac{8}{{3\pi^{2}}}\right)\left(\frac{2}{{d_{1,x}d_{1,y}d_{1,z}}}\left((E-E_{1})^{3}+E_{1}^{3}\right)+\frac{1}{{d_{2,x}d_{2,y}d_{2,z}}}\left((E-E_{2})^{3}+E_{2}^{3}\right)\right).$
(6)

It can be observed from the extracted figures that

$\displaystyle d_{1,x}$ $\displaystyle=3.3\text{ \AA eV},d_{1,y}=1.9\text{ \AA
eV},d_{1,z}=0.3\text{ \AA eV},$ $\displaystyle d_{2,x}$
$\displaystyle=3.0\text{ \AA eV},d_{2,y}=1.9\text{ \AA eV},d_{2,z}=0.2\text{
\AA eV}.$

Therefore, the doping electron concentration in cm-3 of TaP for the desired
Fermi level shift, E, in meV is

$\displaystyle\left(\frac{8\times
10^{15}}{3\pi^{2}}\right)\left(\frac{2}{1.9}\left((E-24)^{3}+24^{3}\right)+\frac{1}{1.1}\left((E+40)^{3}-40^{3}\right)\right).$
(7)

The plot depicting the equation (7) for $E\in[-50,50]$ meV is shown in Fig.
S1c. By substituting $E$ with $E_{1}$, and $E_{2}$, the model predicts the
doping electron concentration of 5.3$\times$1019 cm-3 for W1, and doping hole
concentration of 8.6$\times$1019 cm-3 for W2, respectively.

###  A.3 Ab initio Calculation of H-doped TaP

ab

cd

Figure S2: Convergence test on simulation’s hyperparameters: a kinetic energy
cutoff for wavefunction (ecutwfc), b kinetic energy cutoff for charge density
(ecutrho), c Gaussian spreading for Brillouin-zone integration in metals
(degauss), and d number of k-point grids on each axis.

To estimate the correlation between hydrogen doping concentration in TaP and
the shift in the Fermi level, ab initio calculations were done. It involves
determining the count of unoccupied states situated between the initial and
adjusted Fermi levels. In this work, Quantum Espresso version 7.2 was used as
the ab initio calculation program. The input material structure is obtained
from the MaterialsProject in the primitive cell format [51]. The
pseudopotential files for Ta and P (both with scalar and full relativistic)
are from the Quantum Espresso database. All pseudo potentials are Projector-
Augmented Wave (PAW) types with the Perdew-Burke-Ernzerhof approximation
method for solid (PBEsol). It was initiated by performing convergence test on
simulation hyperparameters: kinetic energy cutoff for wavefunction (ecutwfc),
and charge density (ecutrho), Gaussian spreading for Brillouin-zone
integration in metals (degauss), and number of k-point grid on each axis (nk1,
nk2, nk3). Fig. S2 show the convergence test results.

abc Figure S3: a The electronic band structure of TaP. b The energy dependence
of the density of states (Top), and number of states per unit cell (Bottom,
left-hand axis) that needed to be filled with the doped electron to effectuate
the shift of the Fermi level. The variation of doping electron concentration
is shown by right-hand axis. c A magnified view of b for the range
$E\in[-50,50]$ meV with respect to Fermi level.

Subsequently, structure relaxation was performed on the structure from the
MaterialsProject to get a stable structure according to the pseudopotentials
used. This stable structure was used in Quantum Espresso to calculate the
self-consistent function, and non-self-consistent calculation with high
k-point density grids [52]. Then, we perform density of states and band
structure calculation pathway. The band structure and density of the states
plot are shown in Fig. S3a and Fig. S3b, respectively. The doping relation was
achieved by integrating the density of states starting from the Fermi level.
This allowed the determination of the number of states per unit cell that
needed to be filled with the doped electron to effectuate the shift of the
Fermi level (shown in the bottom of Fig. S3b). Finally, the doping electron
concentration was calculated by dividing the unit cell volume. Our calculation
predicts the doping electron concentration of 9.7$\times$1019 cm-3 for W1, and
doping hole concentration of 1.4$\times$1020 cm-3 for W2.

###  A.4 SRIM and TRIM Calculation

We employed the SRIM (stopping and range of ions in matter) and TRIM
(transport of ions in matter) software [41] to simulate realistic irradiation
conditions and investigate the defect profile. TRIM simulations provided us
with the flexibility to adjust parameters such as target type, thickness, and
incident ion beam energy. Through the TRIM simulations, we obtained valuable
results that revealed the depth profile generated by both the incident ions
and recoils. For our study, we deliberately chose H- ions with an energy of 20
keV to irradiate the TaP crystal. This selection was driven by our specific
interest in doping the ions near the surface of the sample ($\sim$ 5000 Å). By
considering the beam current and aperture size, we calculated the incident ion
flux. Furthermore, we accounted for the subtraction of secondary electrons and
the total duration of irradiation, enabling us to determine the number of
implanted ions in each sample. In our experimental setup, we utilized a 4 mm
diameter beam to irradiate samples of size 2$\times$2 mm2. The beam current
for S3m was 140 nA, and for S20m and S2h was 170 nA. We varied the irradiation
time accordingly to achieve the targeted doping concentration of $\simeq$ 1019
cm-3, 1020 cm-3 and 1021 cm-3 for S3m, S20m and S2h samples as per our
mathematical and computational calculations. The total doping concentration in
the sample is determined by multiplying the irradiation dose for each sample
by the penetration depth. The irradiation dose is summarized in the main text
(Table 1).

###  A.5 Carrier concentration and mobility

abcS3mS20mS2h

defS3mS20mS2h

Figure S4: Magnetic field dependent of longitudinal resistivity, $\rho_{xx}$,
at different temperatures up to 300 K for a S3m, b S20m, and c S2h. d-f The
same as a-c, respectively, but for $\rho_{xy}$.

The electrical measurements of the samples were carried out using the
electrical transport option (ETO) of the physical property measurement system
(PPMS). We utilized the ETO with a six-probe configuration to obtain the
required data, enabling simultaneous measurements of both $\rho_{xx}$ and
$\rho_{xy}$. We employed a symmetric probe configuration to minimize the
influence of manually made contact misalignment. In this setup, the
longitudinal resistivity $\rho_{xx}$ exhibited symmetry with respect to the
applied magnetic field, while the transverse resistivity $\rho_{xy}$ displayed
antisymmetry. We employed the following equations to achieve the symmetric
probe:

$\rho_{xx}(B)=\frac{\rho_{xx}(+B)+\rho_{xx}(-B)}{2},\,\rho_{xy}(B)=\frac{\rho_{xy}(+B)-\rho_{xy}(-B)}{2}.$
(8)

The magnetic field dependence of $\rho_{xx}$ and $\rho_{xy}$ at different
temperatures is illustrated in Fig. S4. Shubnikov-de Haas (SdH) oscillations
can be observed in both $\rho_{xx}$ and $\rho_{xy}$ at low temperatures, up to
15 K for S3m, and 10 K for S20m. However, no such oscillations are observed
for S2h sample, indicating that the presence of high disorder suppresses the
quantum oscillations. Thus, our further experimental analysis primarily
focuses on the S3m, and S20m samples as our main interest lies in precise
Fermi level engineering in a topological Weyl semimetal.

abS3m

cd

ef

Figure S5: Temperature variation of a longitudinal conductivity,
$\sigma_{xx}$, and b transverse conductivity, $\sigma_{xy}$, for S3m. The
fitting results of the curves in a, and b at 2 K with two-band model are shown
in c, and d, respectively. The fitting RMSRE (Root Mean Squared Relative
Error) values at different temperatures for e $\sigma_{xx}$, and f
$\sigma_{xy}$.

abS20m

cd

ef

Figure S6: Temperature variation of a longitudinal conductivity,
$\sigma_{xx}$, and b transverse conductivity, $\sigma_{xy}$, for S20m. The
fitting results of the curves in a, and b at 2 K with two-band model are shown
in c, and d, respectively. The fitting RMSRE (Root Mean Squared Relative
Error) values at different temperatures for e $\sigma_{xx}$, and f
$\sigma_{xy}$.

To investigate the impact of irradiation on carrier concentration and
mobility, we calculated the longitudinal conductivity, $\sigma_{xx}$, and
transverse conductivity, $\sigma_{xy}$, using the following equations:

$\sigma_{xx}=\frac{\rho_{xx}}{{\rho_{xx}^{2}+\rho_{xy}^{2}}},\,\sigma_{xy}=-\frac{\rho_{xy}}{{\rho_{xx}^{2}+\rho_{xy}^{2}}}.$
(9)

The field dependence of $\sigma_{xx}$ and $\sigma_{xy}$ at different
temperatures are illustrated in Fig. S5a and Fig. S5b for S3m, respectively,
and Fig. S6a and Fig. S6b for S20m, respectively. To extract the carrier
concentration and mobility, we simultaneously fit the $\sigma_{xx}$ and
$\sigma_{xy}$ data using two-band model defined by:

$\sigma_{xx}=\frac{n_{e}\mu_{e}e}{1+(\mu_{e}B)^{2}}+\frac{n_{h}\mu_{h}e}{1+(\mu_{h}B)^{2}},\,\sigma_{xy}=\left(\frac{n_{h}\mu_{h}^{2}}{1+(\mu_{h}B)^{2}}-\frac{n_{e}\mu_{e}^{2}}{1+(\mu_{e}B)^{2}}\right)eB$
(10)

Fig. S5c-d, and Fig. S6c-d show example fitting results of the curves for 2 K
temperature for S3m, and S20m, respectively. The root mean squared relative
error (RMSRE) linked to the fitting procedure at different temperatures using
the two-band model is displayed in Fig. S5e-f, and Fig. S6e-f for S3m, and
S20m, respectively. The main text includes the temperature variation of
electron and hole concentration, along with their corresponding mobility, for
each doping level and compares them with pristine TaP (Fig. 4). Our findings
reveal that H--ion implantation significantly impacts the Fermi surface by
considerably elevating the electron-hole compensation temperature despite the
mobility being reduced by an order of magnitude.

###  A.6 Analysis of Quantum Oscillation

Magnetoresistance (MR) data was extracted from $\rho_{xx}$ data, using the
following relation:

$\text{MR}=\frac{\rho_{xx}(B)-\rho_{xx}(0)}{\rho_{xx}(0)}\times 100\%$ (11)

The quantum oscillation data, $\Delta$MR, was obtained by subtracting
background signals from the MR data. To ensure robustness, we employed three
independent methods to extract and analyze $\Delta$MR: (a) background-free
curvature (Fig. S7 and Fig. S8), (b) polynomial background (partially in Fig.
3, and fully in Fig. S9 and Fig. S10), and (c) piece-wise polynomial
background (Fig. S11 and Fig. S12). In (a), a second-order derivative with
respect to the magnetic field $B$ is applied to the MR data, automatically
eliminating all constant, linear, and quadratic terms, which mostly contribute
to background instead of oscillation, without the need for actual background
selection.

Notably, the fast Fourier transform (FFT) analysis unveiled distinct
oscillation frequencies. After applying a standard signal filtering process
through FFT against $1/B$, we isolated the oscillation components from the
data and determined the corresponding Landau levels (LLs). The LL index fan
corresponding to the different frequencies for each fitting model has been
plotted. All three methods consistently support the existence of a low-
frequency carrier pocket, Fα 1.87 T - 3.94 T for S3m and 0.35 T - 1.14 T for
S20m. In contrast, the pristine TaP exhibits the same behavior in the range Fα
2.3 T - 4 T and Fβ 18 T - 19 T [24]. Our comprehensive analysis confirms the
presence of new oscillation frequencies, suggesting that the carrier pockets
of the Weyl node can indeed reach the desired 0 Landau level, and additional
Weyl nodes appear with H- ion irradiation.

abcS3m

def

ghi

Figure S7: MR analysis of S3m with background-free curvature method. MR plots
and their corresponding second derivatives which work as $\Delta$MR at a 2 K,
b 5 K, and c 15 K. d-f Fast Fourier transform (FFT) of a-c, respectively. g-i
display the Shubnikov-de Haas (SdH) oscillation and Landau level fan plot of
d-f, respectively.

abcS20m

def

ghi

Figure S8: MR analysis of S20m with background-free curvature method. MR plots
and their corresponding second derivatives, which work as $\Delta$MR at a 2 K,
b 5 K, and c 10 K. d-f Fast Fourier transform (FFT) of a-c, respectively. g-i
display the Shubnikov-de Haas (SdH) oscillation and Landau level fan plot of
d-f, respectively.

abcS3m

def

ghi

Figure S9: MR analysis of S3m with polynomial background method. MR plots and
$\Delta$MR at a 2 K, b 5 K, and c 15 K. d-f Fast Fourier transform (FFT) of
a-c, respectively. g-i display the Shubnikov-de Haas (SdH) oscillation and
Landau level fan plot of d-f, respectively.

abcS20m

def

ghi

Figure S10: MR analysis of S20m with polynomial background method. MR plots
and $\Delta$MR at a 2 K, b 5 K, and c 10 K. d-f Fast Fourier transform (FFT)
of a-c, respectively. g-i display the Shubnikov-de Haas (SdH) oscillation and
Landau level fan plot of d-f, respectively.

abcS3m

def

ghi

Figure S11: MR analysis of S3m with piece-wise polynomial background method.
MR plots and $\Delta$MR at a 2 K, b 5 K, and c 15 K. d-f Fast Fourier
transform (FFT) of a-c, respectively. g-i display the Shubnikov-de Haas (SdH)
oscillation and Landau level fan plot of d-f, respectively.

abcS20m

def

ghi

Figure S12: MR analysis of S20m with piece-wise polynomial background method.
MR plots and $\Delta$MR at a 2 K, b 5 K, and c 10 K. d-f Fast Fourier
transform (FFT) of a-c, respectively. g-i display the Shubnikov-de Haas (SdH)
oscillation and Landau level fan plot of d-f, respectively.

###  A.7 Extended Kohler’s rule of magnetoresistance

The MR curves do not coincide onto a single curve when plotted against
$B/\rho_{xx}(0)$, indicating that the $\text{MR}=f(B/\rho_{xx}(0))$ rule
proposed by Kohler is not applicable to S3m and S20m, as observed in our
study. However, an intriguing observation is that all the curves in Fig.
S13a-b are nearly parallel with each other. This suggests the possibility of a
single temperature-dependent multiplier affecting MR (y-axis) or
$B/\rho_{xx}(0)$ (x-axis), causing them to coincide or merge into one curve.
To understand this behavior further, MR curves are normalized by a scaling
procedure using the 300 K curve, and adjust the parameter $n_{T}$ individually
for each curve [53]. The outcome is clearly visible in Fig. S13c-d, where all
the curves seamlessly coincide onto the 300 K curve when the data is scaled
using $B/\rho_{xx}(0)n_{T}$. The corresponding values of $n_{T}$ for the
different curves at various temperatures can be found in Fig. S14. This
scaling procedure sheds light on the temperature-dependent behavior of the MR
and provides valuable insights into the underlying physical mechanisms
governing the observed phenomena. We have further normalized $n_{T}$ and
compared it with the data extracted by using a two-band model from the
relation [53]:

$n_{T}=e\left[\frac{\left(\sum_{i}n_{i}\mu_{i}\right)^{3/2}}{\left(\sum_{i}n_{i}\mu_{i}^{3}\right)^{1/2}}\right].$
(12)

abS3mS20m

cdS3mS20m

Figure S13: a-b Kohler’s rule plot of the data for S3m, and S20m,
respectively. c-d Extended Kohler’s rule plot of the data in a-b,
respectively.

abS3mS20m

cdS3mS20m

efS3mS20m

Figure S14: The temperature dependence of the adjusted parameter $n_{T}$
extracted from Extended Kohler’s rule for a S3m, and b S20m. c-d The same as
a-b, but $n_{T}$ are calculated with fitted parameters from the two-band model
using the relation given in equation (12). e-f depict the normalized $n_{T}$
values, enabling a direct comparison between simple extended Kohler’s rule and
fitted two-band model.


remains usable in practice. More precisely, the adjustments require each DP to
perform a public-key encryption and non-interactive zero-knowledge proof to
record an observation, whereas for set union, they simply randomize a
plaintext counter. For communication cost, using the Fiat-Shamir heuristic and
assuming group elements and plaintext elements of $\mathbb{Z}_{q}$ are roughly
equal size, the total cost during input submission should increase from one
plaintext counter (for each counter held by each DP) to $4$ (two each for the
ciphertext and proof) with other phases of the protocol unchanged.

##  9\. Conclusion and Future Work

We present the PSC protocols for securely computing the cardinality of set
union and set intersection across multiple parties. PSC is secure against an
active adversary controlling all but one honest party, provides adaptive
security for the Data Parties, makes the adversary accountable for disrupting
the protocol, and produces differentially-private outputs. Our implementation
of PSC is available for download at https://github.com/GUSecLab/psc. Future
work includes developing the sketch-based approach (Choi et al., 2020) by
modifying generic MPC protocols to provide all our security properties,
applying efficiency improvements in zero-knowledge proofs (Bayer and Groth,
2012; Bünz et al., 2018), and investigating batching broadcasts (Baum et al.,
2016) to improve performance while mainting our termination guarantees.

###### Acknowledgements.

This work has been supported by the National Science Foundation under grant
CNS-1718498, the Canada Research Chairs program, the Office of Naval Research,
and the Callahan Family Professor endowment. The views expressed in this work
are those of the authors and do not necessarily reflect those of the funding
agencies or organizations.

## References

  * (1)
  * Asharov and Lindell (2017) Gilad Asharov and Yehuda Lindell. 2017.  A full proof of the BGW protocol for perfectly secure multiparty computation.  _Journal of Cryptology_ 30, 1 (2017).
  * Baum et al. (2016) Carsten Baum, Emmanuela Orsini, and Peter Scholl. 2016\.  Efficient secure multiparty computation with identifiable abort. In _Theory of Cryptography Conference (TCC)_.
  * Baum et al. (2020) Carsten Baum, Emmanuela Orsini, Peter Scholl, and Eduardo Soria-Vazquez. 2020.  Efficient Constant-Round MPC with Identifiable Abort and Public Verifiability. In _Annual International Cryptology Conference (Crypto)_.
  * Bayer and Groth (2012) Stephanie Bayer and Jens Groth. 2012.  Efficient Zero-knowledge Argument for Correctness of a Shuffle. In _Annual International Conference on the Theory and Applications of Cryptographic Techniques (Eurocrypt)_.
  * Bellare and Rogaway (1993) Mihir Bellare and Phillip Rogaway. 1993.  Random oracles are practical: A paradigm for designing efficient protocols. In _ACM Conference on Computer and Communications Security (CCS)_.
  * Benaloh (1994) Josh Benaloh. 1994\.  Dense probabilistic encryption. In _Workshop on selected areas of cryptography (SAC)_.
  * Blum et al. (2013) Avrim Blum, Katrina Ligett, and Aaron Roth. 2013\.  A learning theory approach to noninteractive database privacy.  _Journal of the ACM (JACM)_ 60, 2 (2013).
  * Brandt (2005) Felix Brandt. 2005\.  Efficient cryptographic protocol design based on distributed El Gamal encryption. In _International Conference on Information Security and Cryptology (ICISC)_.
  * Bünz et al. (2018) Benedikt Bünz, Jonathan Bootle, Dan Boneh, Andrew Poelstra, Pieter Wuille, and Greg Maxwell. 2018.  Bulletproofs: Short proofs for confidential transactions and more. In _IEEE Symposium on Security and Privacy (S &P)_.
  * Camenisch et al. (2017) Jan Camenisch, Manu Drijvers, and Maria Dubovitskaya. 2017\.  Practical UC-secure delegatable credentials with attributes and their application to blockchain. In _ACM Conference on Computer and Communications Security (CCS)_.
  * Canetti (2001) Ran Canetti. 2001\.  Universally Composable Security: A New Paradigm for Cryptographic Protocols. In _Foundations of Computer Science (FOCS)_.
  * Canetti et al. (2002) Ran Canetti, Yehuda Lindell, Rafail Ostrovsky, and Amit Sahai. 2002.  Universally composable two-party and multi-party secure computation. In _Symposium on Theory of Computing (STOC)_.
  * Canetti and Rabin (2002) Ran Canetti and Tal Rabin. 2002.  Universal Composition with Joint State.  Cryptology ePrint Archive, Report 2002/047.  https://eprint.iacr.org/2002/047.
  * Choi et al. (2020) Seung Geol Choi, Dana Dachman-Soled, Mukul Kulkarni, and Arkady Yerukhimovich. 2020.  Differentially-Private Multi-Party Sketching for Large-Scale Statistics.  _Proceedings on Privacy Enhancing Technologies_ 3 (2020).
  * Damgård (2010) Ivan Damgård. 2010\.  On $\Sigma$-protocols.  Lecture Notes on Cryptologic Protocol Theory, v.2.
  * Damgård et al. (2012) Ivan Damgård, Valerio Pastro, Nigel Smart, and Sarah Zakarias. 2012.  Multiparty computation from somewhat homomorphic encryption. In _Annual International Cryptology Conference (Crypto)_.
  * De Cristofaro et al. (2012) Emiliano De Cristofaro, Paolo Gasti, and Gene Tsudik. 2012\.  Fast and private computation of cardinality of set intersection and union. In _International Conference on Cryptology and Network Security (CANS)_.
  * de Kok (2020) Daniël de Kok. 2020\.  Go par package for parallel for-loops.  https://github.com/danieldk/par.
  * Dingledine et al. (2004) Roger Dingledine, Nick Mathewson, and Paul Syverson. 2004\.  Tor: The Second-Generation Onion Router. In _USENIX Security Symposium (USENIX)_.
  * Dolev and Strong (1983) Danny Dolev and H. Raymond Strong. 1983.  Authenticated Algorithms for Byzantine Agreement.  _SIAM J. Comput._ 12, 4 (1983).
  * Durand and Flajolet (2003) Marianne Durand and Philippe Flajolet. 2003.  Loglog counting of large cardinalities. In _European Symposium on Algorithms (ESA)_.
  * Dwork et al. (2006a) Cynthia Dwork, Krishnaram Kenthapadi, Frank McSherry, Ilya Mironov, and Moni Naor. 2006a.  Our Data, Ourselves: Privacy via Distributed Noise Generation. In _Advances in Cryptology (Eurocrypt)_.
  * Dwork et al. (2006b) Cynthia Dwork, Frank McSherry, Kobbi Nissim, and Adam Smith. 2006b.  Calibrating Noise to Sensitivity in Private Data Analysis. In _Theory of Cryptography Conference (TCC)_.
  * Dwork et al. (2014) Cynthia Dwork, Aaron Roth, et al. 2014\.  The algorithmic foundations of differential privacy.  _Foundations and Trends in Theoretical Computer Science_ 9, 3-4 (2014).
  * Egert et al. (2015) Rolf Egert, Marc Fischlin, David Gens, Sven Jacob, Matthias Senker, and Jörn Tillmanns. 2015.  Privately Computing Set-Union and Set-Intersection Cardinality via Bloom Filters. In _Australasian Conference on Information Security and Privacy_.
  * Elahi et al. (2014) Tariq Elahi, George Danezis, and Ian Goldberg. 2014\.  PrivEx: Private Collection of Traffic Statistics for Anonymous Communication Networks. In _ACM Conference on Computer and Communications Security (CCS)_.
  * Fenske et al. (2017) Ellis Fenske, Akshaya Mani, Aaron Johnson, and Micah Sherr. 2017\.  Distributed Measurement with Private Set-Union Cardinality. In _ACM Conference on Computer and Communications Security (CCS)_. ACM.
  * Fiat and Shamir (1987) Amos Fiat and Adi Shamir. 1987.  How to Prove Yourself: Practical Solutions to Identification and Signature Problems. In _Advances in Cryptology (CRYPTO ’86)_.
  * Freedman et al. (2004) Michael J Freedman, Kobbi Nissim, and Benny Pinkas. 2004\.  Efficient Private Matching and Set Intersection. In _Advances in Cryptology (Eurocrypt)_.
  * Furukawa et al. (2003) Jun Furukawa, Hiroshi Miyauchi, Kengo Mori, Satoshi Obana, and Kazue Sako. 2003.  An Implementation of a Universally Verifiable Electronic Voting Scheme Based on Shuffling. In _Financial Cryptography_ _(FC’02)_.
  * Goldreich (2001) O. Goldreich. 2001\.  _Foundations of Cryptography: Volume 2, Basic Applications_.  Cambridge University Press.
  * Goldreich et al. (1987) Oded Goldreich, Silvio Micali, and Avi Wigderson. 1987\.  How to play ANY mental game. In _ACM Symposium on Theory of Computing (STOC)_.
  * Goldreich and Oren (1994) Oded Goldreich and Yair Oren. 1994.  Definitions and properties of zero-knowledge proof systems.  _Journal of Cryptology_ 7, 1 (1994).
  * Groth (2003) Jens Groth. 2003\.  A Verifiable Secret Shuffle of Homomorphic Encryptions. In _Theory and Practice in Public Key Cryptography (PKC)_.
  * Groth et al. (2006) Jens Groth, Rafail Ostrovsky, and Amit Sahai. 2006\.  Perfect non-interactive zero knowledge for NP. In _Annual International Conference on the Theory and Applications of Cryptographic Techniques (EUROCRYPT)_.
  * Hazay et al. (2012) Carmit Hazay, Gert Læssøe Mikkelsen, Tal Rabin, Tomas Toft, and Angelo Agatino Nicolosi. 2012.  Efficient RSA Key Generation and Threshold Paillier in the Two-Party Setting. In _Topics in Cryptology – CT-RSA_.
  * Hazay and Nissim (2012) Carmit Hazay and Kobbi Nissim. 2012.  Efficient Set Operations in the Presence of Malicious Adversaries.  _Journal of Cryptology_ 25, 3 (2012).
  * Inan et al. (2010) Ali Inan, Murat Kantarcioglu, Gabriel Ghinita, and Elisa Bertino. 2010.  Private record matching using differential privacy. In _International Conference on Extending Database Technology_.
  * Ishai et al. (2014) Yuval Ishai, Rafail Ostrovsky, and Vassilis Zikas. 2014\.  Secure multi-party computation with identifiable abort. In _Annual Cryptology Conference (CRYPTO)_.
  * Jansen and Johnson (2016) Rob Jansen and Aaron Johnson. 2016.  Safely Measuring Tor. In _ACM Conference on Computer and Communications Security (CCS)_.
  * Kasiviswanathan and Smith (2014) Shiva P Kasiviswanathan and Adam Smith. 2014.  On the ‘Semantics’ of Differential Privacy: A Bayesian Formulation.  _Journal of Privacy and Confidentiality_ 6, 1 (2014).
  * Katz et al. (2013) Jonathan Katz, Ueli Maurer, Björn Tackmann, and Vassilis Zikas. 2013.  Universally composable synchronous computation. In _Theory of Cryptography Conference (TCC)_.
  * Kiayias et al. (2016) Aggelos Kiayias, Hong-Sheng Zhou, and Vassilis Zikas. 2016\.  Fair and Robust Multi-party Computation Using a Global Transaction Ledger. In _Advances in Cryptology (EUROCRYPT)_.
  * Kissner and Song (2005) Lea Kissner and Dawn Song. 2005.  Privacy-Preserving Set Operations. In _Annual International Cryptology Conference (Crypto)_.
  * kyber. (2020) kyber. 2020.  kyber: DEDIS Advanced Crypto Library for Go.  https://godoc.org/go.dedis.ch/kyber.
  * Larraia et al. (2014) Enrique Larraia, Emmanuela Orsini, and Nigel P Smart. 2014\.  Dishonest Majority Multi-Party Computation for Binary Circuits. In _Annual International Cryptology Conference (Crypto)_.
  * Lindell (2005) Yehida Lindell. 2005\.  Secure multiparty computation for privacy preserving data mining.  In _Encyclopedia of Data Warehousing and Mining_. 1005–1009.
  * Lindell (2015) Yehuda Lindell. 2015\.  An Efficient Transform from Sigma Protocols to NIZK with a CRS and Non-programmable Random Oracle. In _Theory of Cryptography (TCC)_.
  * Lindell et al. (2015) Yehuda Lindell, Benny Pinkas, Nigel P Smart, and Avishay Yanai. 2015.  Efficient constant round multi-party computation combining BMR and SPDZ. In _Annual Cryptology Conference (Crypto)_.
  * Mani and Sherr (2017) Akshaya Mani and M. Sherr. 2017.  Histor$\epsilon$: Differentially Private and Robust Statistics Collection for Tor. In _Network and Distributed System Security Symposium (NDSS)_.
  * McCoy et al. (2008) Damon McCoy, Kevin Bauer, Dirk Grunwald, Tadayoshi Kohno, and Douglas Sicker. 2008.  Shining Light in Dark Places: Understanding the Tor Network. In _Privacy Enhancing Technologies Symposium (PETS)_.
  * McSherry and Talwar (2007) Frank McSherry and Kunal Talwar. 2007.  Mechanism design via differential privacy. In _Foundations of Computer Science (FOCS)_.
  * Melis et al. (2016) Luca Melis, George Danezis, and Emiliano De Cristofaro. 2016\.  Efficient Private Statistics with Succinct Sketches. In _Network and Distributed System Security Symposium (NDSS)_.
  * Neff (2001) C. Andrew Neff. 2001\.  A Verifiable Secret Shuffle and its Application to E-voting. In _ACM Conference on Computer and Communications Security (CCS)_.
  * Nguyen et al. (2004) Lan Nguyen, Rei Safavi-Naini, and Kaoru Kurosawa. 2004\.  Verifiable Shuffles: A Formal Model and a Paillier-based Efficient Construction with Provable Security. In _Applied Cryptography and Network Security (ACNS)_.
  * Nissim et al. (2007) Kobbi Nissim, Sofya Raskhodnikova, and Adam Smith. 2007\.  Smooth Sensitivity and Sampling in Private Data Analysis. In _Symposium on Theory of Computing (STOC)_.
  * Partridge and Allman (2016) Craig Partridge and Mark Allman. 2016.  Ethical considerations in network measurement papers.  _Commun. ACM_ 59, 10 (2016).
  * Pettai and Laud (2015) Martin Pettai and Peeter Laud. 2015.  Combining differential privacy and secure multiparty computation. In _Annual Computer Security Applications Conference (ACSAC)_.
  * Ran Canetti (2004) Shai Halevi Ran Canetti, Oded Goldreich. 2004\.  The Random Oracle Methodology, Revisited.  _Journal of the ACM (JACM)_ 51, 4 (2004).
  * Schnorr (1991) Claus-Peter Schnorr. 1991\.  Efficient Signature Generation by Smart Cards.  _Journal of Cryptology_ 4, 3 (1991).
  * Soghoian (2011) Christopher Soghoian. 2011\.  Enforced Community Standards For Research on Users of the Tor Anonymity Network. In _Workshop on Ethics in Computer Security Research (WECSR)_.
  * Stanojevic et al. (2017) Rade Stanojevic, Mohamed Nabeel, and Ting Yu. 2017.  Distributed Cardinality Estimation of Set Operations with Differential Privacy. In _IEEE Symposium on Privacy-Aware Computing (PAC)_.
  * Stoica et al. (2001) Ion Stoica, Robert Morris, David Karger, M Frans Kaashoek, and Hari Balakrishnan. 2001.  Chord: A Scalable Peer-to-Peer Lookup Service for Internet Applications. In _Conference on Applications, Technologies, Architectures, and Protocols for Computer Communications (SIGCOMM)_.
  * Struik (2021) Rene Struik. 2021\.  _Alternative Elliptic Curve Representations_.  Internet-Draft draft-ietf-lwig-curve-representations-20. Internet Engineering Task Force.  Work in Progress.
  * Tor Metrics. (2020) Tor Metrics. 2020.  https://metrics.torproject.org/.
  * Vaidya and Clifton (2005) Jaideep Vaidya and Chris Clifton. 2005.  Secure Set Intersection Cardinality with Application to Association Rule Mining.  _Journal of Computer Security_ 13, 4 (2005).
  * Varda (2008) Kenton Varda. 2008\.  Protocol Buffers: Google’s Data Interchange Format.  https://opensource.googleblog.com/2008/07/protocol-buffers-googles-data.html.
  * Wails et al. (2019) Ryan Wails, Aaron Johnson, Daniel Starin, Arkady Yerukhimovich, and S. Dov Gordon. 2019\.  Stormy: Statistics in Tor by Measuring Securely.  _ACM Conference on Computer and Communications Security (CCS)_.
  * Wang et al. (2017) Xiao Wang, Samuel Ranellucci, and Jonathan Katz. 2017\.  Global-scale secure multiparty computation. In _ACM Conference on Computer and Communications Security (CCS)_.
  * Wikström (2005) Douglas Wikström. 2005\.  A Sender Verifiable Mix-Net and a New Proof of a Shuffle. In _International Conference on the Theory and Application of Cryptology and Information Security (ASIACRYPT)_.

##  Appendix A Commitments and Coins

For reference, in the section we repeat the one-to-many commitment
functionality given in (Canetti et al., 2002) (Figure 12), present an
accountable coin-flipping functionality (Figure 13), and give a protocol that
UC-realizes the coin flipping functionality and the corresponding proof
(Figure 14).

$\mathcal{F}_{\textsc{MCOM}}^{\textsc{1:M}}$

  1. (1)

Commit Phase: Upon receiving a message (commit, $m$) from $\textrm{CP}_{j}$
with $m\in\\{0,1\\}$, record the tuple ($\textrm{CP}_{j}$, $m$) and send
(receipt, $\textrm{CP}_{j}$) to all parties.

  2. (2)

Prove Phase: Upon receiving a message (reveal) from $\textrm{CP}_{j}$, if a
tuple ($\textrm{CP}_{j}$, $m$) was previously recorded, then send the message
(reveal, $m$) to all parties.

Figure 12. $\mathcal{F}_{\textsc{MCOM}}^{\textsc{1:M}}$, the ideal
functionality for one-to-many commitments.

$\mathcal{F}_{\textsc{COIN}}^{\textsc{}}$

  1. (1)

(Round 1) Wait to receive GenerateCoins from each CP

  2. (2)

(Round 2) If received GenerateCoins from each CP, generate
$\vec{b}\in\\{0,1\\}^{n}$ at random and broadcast $\vec{b}$ to all CPs.
Otherwise, broadcast (Blame, $L$) where $L$ where $L$ is the set of the CPs
from whom messages are missing from Round 1.

Figure 13. $\mathcal{F}_{\textsc{COIN}}^{\textsc{}}$, the ideal coin flipping
functionality.

$\pi_{\textsc{COIN}}$

  1. (1)

(Round 1) Commitments.  Each $\textrm{CP}_{j}$ selects $b_{j}\in\\{0,1\\}$ at
random and commits to $b_{j}$ to every other CP using
$\mathcal{F}_{\textsc{MCOM}}^{\textsc{1:M}}$.

  2. (2)

If a CP has received a commitment receipt from all other CPs, each addressed
to the set of all CPs, proceed to Round 2. Otherwise, broadcast (Blame, $L$)
where $L$ is the set of CPs for which a correct receipt was not received and
then exit.

  3. (3)

(Round 2) Opening.  Each $\textrm{CP}_{j}$ opens $b_{j}$ to all CPs using
$\mathcal{F}_{\textsc{MCOM}}^{\textsc{1:M}}$. If a CP receives a valid opening
from each CP, output $b=\bigoplus_{j}b_{j}$. Otherwise, output (Blame, $L$)
where $L$ is the set of CPs who have failed to open their commitments to a
single bit.

Figure 14. $\pi_{\textsc{COIN}}$, an accountable coin-flipping protocol

###### Theorem A.1.

The protocol $\pi_{\textsc{COIN}}$ UC-realizes the
$\mathcal{F}_{\textsc{COIN}}^{\textsc{}}$ functionality in the static-
corruption, synchronous, ($\mathcal{F}_{\textsc{BC}}^{\textsc{}}$,
$\mathcal{F}_{\textsc{MCOM}}^{\textsc{1:M}}$)-hybrid model if there is at
least one honest CP.

###### Proof.

We define a simulator $S$ as follows: $S$ runs a copy of the adversary $A$,
and simulates the honest parties and ideal functionalities honestly. If a
malicious party commits during the commit phase and opens it to an element of
$\\{0,1\\}$ in the opening phase, S sends GenerateCoin to
$\mathcal{F}_{\textsc{COIN}}^{\textsc{}}$ on behalf of that malicious party.
We observe that every dishonest party $P$ either completes a valid opening in
the second phase of the protocol, in which case
$\mathcal{F}_{\textsc{COIN}}^{\textsc{}}$ has received GenerateCoin from $P$
(through $S$), or $P$ has failed to produce a valid opening, in which it is
blamed by all honest (simulated) CPs. Finally, in the case that all parties
open commitments and $\mathcal{F}_{\textsc{COIN}}^{\textsc{}}$ outputs a
uniformly random bit, this bit is identically distributed to $\oplus_{j}b_{j}$
in the simulation, since $b_{j}$ is random for honest $j$, and dishonest CPs
must commit to their bits without any information about the honest parties’
selections. ∎

We will need a vector of $n$ coins, so we write
$\mathcal{F}_{\textsc{N-COIN}}^{\textsc{}}$ as $n$ independent copies of
$\mathcal{F}_{\textsc{COIN}}^{\textsc{}}$, which can be joined into a single
protocol using the JUC theorem from (Canetti and Rabin, 2002). Similarly, we
will require commitments to strings of bits, so we write
$\mathcal{F}_{\textsc{S-MCOM}}^{\textsc{1:M}}$ to indicate the protocol where
$\mathcal{F}_{\textsc{MCOM}}^{\textsc{1:M}}$ is invoked multiple times to
produce a commitment to a bit string rather than a single bit.