# CoMP-Enabled RAN Slicing for Tactile Internet

Peng Yang, Xing Xi, Tony Q. S. Quek, and Hyundong Shin

###### Abstract

Tactile Internet (TI) enables the omnipresence and exchange of tactile
experiences across the globe via the ultra-reliable and ultra-responsive
connectivity. This article argues for coordinated multi-point (CoMP) enabled
radio access network (RAN) slicing as an efficient solution that satisfies the
stringent reliable and responsive connectivity requirements for supporting
tactile applications. This article presents the emerging challenges when
accommodating CoMP-enabled RAN slicing in the TI ecosystem and expounds on the
functional split of a CoMP-enabled RAN. The implementation prototype of CoMP-
enabled RAN slicing for TI is also elaborated when the tactile applications
coexist with other vertical applications. Finally, this article studies a use
case of enabling TI-included application multiplexing as an example of CoMP-
enabled RAN slicing for TI.

##  I Introduction

The Tactile Internet (TI), poised to shift networks from content transmission
to skill-set transmission, is envisioned to revolutionize most segments of our
society [1]. To realize the grand envision, TI is expected to deliver remote
control and physical haptic information in real time via the ultra-reliable
and ultra-responsive (i.e., the end-to-end (E2E) latency on the order of one
millisecond) network connectivity [2]. Owing to the stringent reliability and
response requirements, one might be persuaded to construct a network
separately designed for E2E control and physical haptic communications. This
scheme, however, is infeasible in terms of the capital (CAPEX) and operational
(OPEX) expenditures. To improve the cost efficiency of providing E2E tactile
services, TI should be established under the existing and extending physical
network infrastructure. Furthermore, there is a consensus that existing
networks must evolve in the direction of simultaneously supporting multiple
vertical applications (e.g., haptic communications, massive machine type
communications (mMTC), and enhanced mobile broadband (eMBB) [3]) with
completely different requirements on reliability, latency, throughput, and
availability [1].

###  I-A Definition of CoMP-enabled RAN slicing

Radio access network (RAN) slicing is a crucial concept of enabling the
existing network infrastructure to be efficiently shared among diverse
vertical applications [4]. The concept of RAN slicing is defined for running
multiple logical or virtual networks as independent business operations on a
common RAN. Each network slice indicates an independent logical network
tailored to accommodate the quality-of-service (QoS) requirement of a
particular application. Besides, to achieve the hard connectivity
requirements, some advanced radio access technologies (RATs) such as massive
MIMO, coordinated multi-point (CoMP) transmission, full-duplex are required to
be explored [2]. Especially, CoMP transmission, creating spatial diversity
with redundant communication paths, has been extensively considered as one of
the most significant RATs that fulfill the hard-bound reliability requirement
for tactile applications [1, 5]. Hence, the concept of CoMP-enabled RAN
slicing is defined as the partition of a RAN exploring the CoMP transmission
technique (called CoMP-enabled RAN) into multiple virtual networks.
Nevertheless, adapting CoMP-enabled RAN slicing to the TI ecosystem brings
many potential challenges, in part owing to the difficulties in abstracting
the CoMP-enabled RAN into multiple slices, as described below.

###  I-B Design challenges

Effective network resource allocation: Network slicing will quickly run into
the limited availability of physical network resources (e.g., transmit power
and spectrum). This limitation is aggravated when network slice resources are
allocated based on tactile traffic with fully predictable arrival time and
amount. If the arrival time and the number of tactile packets cannot be
accurately predicted (e.g., in telesurgery) [1], effectively allocating
network resources to guarantee ultra-reliable control and haptic information
transmission will be challenging. Besides, to satisfy the stringent low-
latency requirement, tactile packets must be immediately transported once
generated. Tactile packets will be blocked and even dropped if not enough
resources are allocated for the slice carrying tactile applications (or TI
slice for brevity), in turn reducing the reliability.

Multi-timescale issue: Network slicing operations (e.g., activate and release
network slices) will be performed in the timescale of minutes to hours to keep
in sync with the slicing of upper layers. In upper-layer slicing processes,
some protocols (e.g., RAN protocol stacks) and functions (e.g., radio resource
control (R2C)) will be activated and configured, the operations of which are
time-consuming [4]. Nevertheless, wireless channels may change on the order of
a millisecond to seconds. Remarkably, the time interval during which the
wireless channel remains unchanged may be much shorter than the RAN slicing
operation duration. Then, optimally conducting time-consuming RAN slicing
operations in rapidly varying wireless channels, which is the so-called multi-
timescale issue of network slicing, should be tackled and is challenging.

Coordination latency reduction: CoMP transmission acting as a spatial
diversity technique that has been standardized [5] can significantly improve
transmission reliability under the low-latency constraint [1]. Further, it can
be used to avoid retransmissions through the message (e.g., payload, command,
and channel state information (CSI)) exchange among all remote radio heads
(RRHs) for coordination. However, the message exchange results in non-
negligible latency, especially in scenarios with low-latency requirements. How
to reduce coordination latency in CoMP transmission is a big challenge.

###  I-C Related work and contributions

####  I-C1 Related work

The concept of network slicing has been extensively investigated in the
literature, wherein a dedicated portion of RAN elements is fully reserved for
specific services. Nevertheless, with TI’s advent, network slicing for TI has
evolved in more agile network frameworks, aiming to attain immersive haptic
interactions and enable varying levels of human-in/out-of-the loop services.
For instance, in [6], the problem of placing virtualized network functions
(VNFs) in a VNF chain, which constituted network slices for serving tactile
applications, and routing tactile traffic through the chain to fulfill the
hard low-latency mandate was studied. Based on radio-aware software-defined
networking (SDN) and NFV techniques, the work in [1] developed a flexible
network architecture of creating different network slices to serve TI-included
multiple applications on common programmable physical infrastructure. In [7],
a two-level medium access control (MAC) scheduling framework for reducing
transmission scheduling latency in network slices, which were constructed to
enable TI-involved application multiplexing over a shared RAN, was developed.
Besides, an elastic cloud-based RAN (CRAN) architecture providing multi-
connectivity for tactile users via dynamically splitting network slices was
proposed in [8].

####  I-C2 Our contributions

In this article, we set forth the fundamental pillars for the exploration of
the CoMP-enabled RAN slicing concept in TI in detail based on the 3GPP
recommended centralized unit (CU) and distributed unit (DU) structure for the
next generation RAN architecture [9]. Specific focus is put on the basic
architectural principles for incorporating CoMP-enabled RAN slicing into the
TI ecosystem. In this regard, this article underlines the crucial domains of
TI architecture, putting emphasis on the functional split of a CoMP-enabled
RAN for TI. Besides, this article elaborates on the implementation prototype
of CoMP-enabled RAN slicing that enables the multiplexing of TI-involved
multiple vertical applications.

Figure 1: TI architecture and functional split of a CoMP-enabled RAN.

##  II TI architecture and functional split

TI offers the medium for delivering touch and actuation in real time. The
haptic sense occurs bilaterally; that is, haptic information is sensed by
imposing an action on an unknown environment and feeling the environment by a
reaction force or a distortion. In this way, TI realizes the haptic control
via the Internet, which is quite different from the conventional Internet. For
the conventional Internet, the feedback can only be audio/visual information.
The exploration on the TI architecture is a crucial problem of enabling TI.
Figure 1 depicts an E2E architecture of TI, which is divided into three
different domains; that is, a master domain, a controlled domain, and a
network domain [1].

###  II-A Master domain

The master domain comprises a tactile operator (human) and a human-system
interface (HSI) that can be a haptic device. The HSI is used to convert the
tactile operator’s inputs to haptic inputs via tactile coding techniques.
Through operating (e.g., touch, feel, and manipulate) the HSI, the tactile
operator can control objects in the controlled domain. Except for outputting
haptic information, the master domain has the provisioning of visual/audio
sensory feedback owing to multi-dimensional characteristics of human tactile
perception. By integrating diverse sensory modalities, the human perception
performance can be significantly improved, and thereby the haptic and non-
haptic feedback greatly help the tactile operator realize the immersive
control and steering of the controlled domain.

###  II-B Controlled domain

The controlled domain mainly includes a controlled robot(s) responding to
tactile operator’s commands to control various objects in the unknown
environment coordinately. The robot reports both cutaneous and kinaesthetic
data and leverages codec to report haptic feedback to the tactile operator.
Thus, the energy between the master and controlled domains is exchanged
through command and haptic (and audio/visual) feedback signal delivery.

###  II-C Network domain

The network domain constructs a network architecture for bilateral
communications between the master and controlled domains. It aims to provide
the tactile operator with an immersive experience in the remote environment.
To this aim, TI needs the ultra-reliable and ultra-responsive network
connectivity, as mentioned at the beginning. The network architecture
comprising a CN and a CoMP-enabled RAN and the correspondingly functional
split, as shown in Figure 1, is expected to be designed to fulfill the
critical connectivity requirements.

The CN is in charge of security processing, defining the access levels of
services available to diverse applications, and handling connection and
mobility management tasks. Accordingly, the CN comprises some dedicated
hardware or virtualized modules (e.g., AUSF, AMF, and home subscribe server
(HSS)). Significant functions of the CoMP-enabled RAN comprise ultra-reliable
and low-latency packet delivery, cost-efficient QoS aware scheduling, and
efficient resource management. Hence, the CoMP-enabled RAN is split into a CU,
a DU pool, and many RUs. This paradigm of functional split also leads to cost
reduction and dynamic scalability via virtualization. Particularly, RUs (or
RRHs) comprise radio-frequency circuitry indicating that RUs’ functionalities
are implemented as physical network functions (PNFs). Other functionalities
distributed over the CU and the DU pool can be virtually implemented as VNFs.
The CU domain contains higher-layer functionalities, and the DU pool includes
lower-layer functionalities. For example, high-level MAC scheduling decisions
(i.e., CoMP) are made in the CU while the time-critical MAC processing (e.g.,
HARQ) resides in the DU pool.

There is a consensus that future networks (including the CoMP-enabled RAN)
must be flexibly designed to be efficiently shared by diverse vertical
applications from both academia and industry. The efficient network sharing
can be achieved by network slicing with diverse network slices being assigned
to serve different types of vertical applications. We next present the
implementation prototype of CoMP-enabled RAN slicing for TI with the
coexistence of other vertical applications and study a use case of application
multiplexing.

##  III Implementation prototype and case study

###  III-A Implementation prototype of CoMP-enabled RAN slicing for TI

Figure 2: Schematic of the implementation prototype of CoMP-enabled RAN
slicing for TI.

Figure 2 illustrates how to implement the CoMP-enabled RAN slicing for TI, the
prototype of which comprises three major components [10]: Software-defined
CoMP RAN coordinator (SRAN-C): calculate and update accepted network slices;
Network slice management function (NSMF): a virtualized function used to
manage as well as control network slices; CN control function (C2F): aim at
configuring the CN based on slices’ QoS requirements.

The SRAN-C module comprises five parts: radio resource control and management
(R2CM), DU pool, CU, user plane anchor (UP-Anchor), and algorithm. It is time-
consuming to create and configure network slices. If some slices are created,
then R2CM needs to configure RAN protocol stacks based on network slices’ QoS
requirements. For instance, some schemes like shortened orthogonal frequency
division multiplexing (OFDM) symbol duration, lower frame error rates, and
CoMP transmission are desired to be explored for TI slices. Inevitably, the
exploration of the CoMP transmission technique results in undesired
coordination latency. Yet, we can resort to the scalable numerology (e.g., the
elastic number of OFDM symbols in a transmit time interval (TTI) and
subcarrier spacing) to bridge the coordination latency to some extent.
Additionally, common/dedicated R2C functions in R2CM will be activated for
application-specific control of radio resources based on VNF/PNF [4]. The DU
pool is introduced to implement PHY collaboration facilitating the CoMP
deployment. MAC scheduler aims at scheduling traffic according to network
conditions to relieve network congestion. Depending on particular
applications, CU will configure and tailor corresponding UP protocol stacks,
for instance, for TI slices, Internet protocol (IP) and header compression
techniques will be disabled, and the RLC function will be transparently
enabled in CU [4]. Due to the requirement for high capacity (data rates)
between the CU and the DUs, optical fiber communications are ideal or non-
constrained fronthaul. Some lower-layer protocols like the open base station
architecture initiative (OBSAI) and the common public radio interface (CPRI)
[11] with a goal of unifying the fronthaul interface also need to be adopted.
UP-Anchor aims to distribute traffic on the basis of the configured network
slice policy and conduct encryption with slice-specific keys as well. For
instance, policy requirements for reliability and security need to be
identified on priority for TI slices. Besides, according to network slice
requests and available resources, the algorithm repetitively calculates,
updates, and reconfigures accepted network slices to obtain the maximum slice
utility. The SRAN-C module also realizes the RAN-CN mapping that is done by
leveraging some configuration protocols (e.g., FlexRAN) [12].

The NSMF deployed in physical/virtualized network infrastructure is in charge
of creating, maintaining, activating, configuring, and releasing network
slices during their life cycles. Through dedicated/common network slice
control functions, NSMF generates a network slice blueprint (i.e., a template)
for every network slice. The network slice blueprint includes the structure,
configuration, control signals, and service flows of a network slice, which
are used to instantiate and control the network slice. A network slice
instance comprises a collection of network functions and resources configured
to fulfill an application’s QoS requirement.

C2F will interpret a network slice blueprint when the network slice request
arrives. Specifically, it will configure the network based on the interpreted
blueprint and identify the optimal paths and servers to deploy VNFs such that
the slice’s QoS requirement can be satisfied. It also reserves several data
centers for supporting network slices. Every data center comprises some
servers (e.g., HSS) to provide VNFs’ services (e.g., identity, independent
subscription and session) for a network slice. Data centers are connected
through high-capacity backhaul and can offer services separately or jointly.

###  III-B Case study of CoMP-enabled RAN slicing: TI-involved application
multiplexing

Based on the implementation prototype, we provide a case study of CoMP-enabled
RAN slicing for TI, mMTC, and eMBB multiplexing to show its potentials in
application multiplexing. Recall that the solution to the challenging issues
of effective network resource allocation and multi-timescale are inevitable
when investigating the CoMP-enabled RAN slicing. Therefore, we formulate the
CoMP-enabled RAN slicing for TI, mMTC, and eMBB multiplexing as a multi-
timescale resource optimization problem with a goal of maximizing the overall
slice utilities. For each type of applications, its slice utility is defined
to reflect users’ QoS mandates in the application and RUs’ total transmit
power consumption. Particularly, we use the E2E latency, codeword error
decoding probability, and packet blocking probability to characterize the QoS
requirements of tactile robots [10]. The achievable data rates and random
access (RA) success probabilities of users are used to measure eMBB and mMTC
users’ QoS requirements, respectively [13]. Besides, in this problem, the
system bandwidth allocated to different network slices is optimized every ten
minutes, and the RU transmit power is optimized every ten seconds. This
problem is constrained by i) total system bandwidth; ii) the maximum RU
transmit power; iii) QoS requirements of tactile robots and eMBB and mMTC
users.

The formulated multi-timescale optimization problem is confirmed to be mixed-
integer non-convex. By exploring an alternating direction method of
multipliers (ADMM) [10], we propose a novel iterative resource allocation
(IRA-ADMM) algorithm to solve the problem. Specifically, by exploring sample
average approximations (SAA), IRA-ADMM converts the multi-timescale problem to
multiple mixed-integer non-convex subproblems of single timescale. To make the
subproblems tractable, an iterative optimization scheme is firstly explored to
handle the mixed-integer issue. Next, a semidefinite relaxation (SDR) scheme
[10] is applied to transform the non-convex subproblems into convex ones. The
solutions (i.e., the bandwidth allocated to each network slice and RU transmit
power) to the single timescale subproblems can then be obtained. Finally, IRA-
ADMM aggregates the achieved slice bandwidth in the single timescale into the
multi-timescale one using the ADMM.

Next, we design a simulation to verify the effectiveness of IRA-ADMM and
compare IRA-ADMM with an iterative resource allocation (IRA) algorithm.
Differing from IRA-ADMM, IRA does not recover the multi-timescale slice
bandwidth from obtained single-timescale slice bandwidth results. We consider
a one km2 communication area wherein all tactile robots and eMBB and mMTC
users are randomly and uniformly distributed, and three RUs are uniformly
deployed to serve these robots and users coordinately. There are two TI
slices, three eMBB slices, and three mMTC slices. The first TI slice includes
three tactile robots with a latency requirement of one millisecond. The second
TI slice contains five robots with a latency requirement of two milliseconds.
The codeword error decoding probability and packet blocking probability
requirements of all tactile robots are set to be $10^{-5}$ and $2\times
10^{-8}$, respectively. Tactile packets are generated by a Poisson
distribution [14]. The number of users in three eMBB slices is four, six, and
eight with the achievable data rate requirements of 6 Mb/s, 4 Mb/s, and 2
Mb/s, respectively. There are $600$ mMTC users in each mMTC slice with each
user’s RA success probability requirement being $0.5$. Other parameter
settings are as follows: noise power is $-110$ dBm, total system bandwidth is
$4$ MHz, the maximum transmit power of each RU is $1$ W, each RU is equipped
with two antennas, all robots and users are single antenna-equipped, the
length of a tactile packet is $160$ bits [10, 13].

Figure 3: Bandwidth allocated to TI slices vs. intensity of arrival tactile
packets.

Figure 3 shows the relationship between the bandwidth allocated to TI slices
and the tactile packet arrival rate when considering the TI and eMBB
multiplexing. We observe that both IRA-ADMM and IRA will allocate the total
system bandwidth to TI slices to guarantee the QoS performance of TI slices
when terminating eMBB services. When tactile and eMBB applications share the
network infrastructure, the network bandwidth allocated to TI slices
monotonously increases with an increasing tactile packet arrival rate. This
observation indicates that the proposed algorithm can effectively solve the
resource allocation issue.

Figure 4: Total TI and eMBB slice utilities vs. intensity of arrival tactile
packets.

Figure 4 plots the trend of the total TI and eMBB slice utilities over
different tactile packet arrival rates. We observe that the obtained total TI
and eMBB slice utilities of IRA-ADMM and IRA monotonously decrease with an
increasing tactile packet arrival rate. Besides, IRA-ADMM achieves greater
total slice utilities than IRA. As an increasing tactile packet arrival rate
results in a decreasing amount of network bandwidth allocated to eMBB slices,
greater RU transmit power should be configured to ensure QoS requirements of
eMBB slices. Thus, a small eMBB slice utility incurring the reduction of the
total slice utilities is obtained. IRA-ADMM obtains greater TI slice utility
than IRA because IRA-ADMM suggests to allocate more bandwidth to TI slices, as
shown in Figure 3.

Figure 5: Total TI and mMTC slice utilities vs. intensity of arrival tactile
packets.

Figure 5 illustrates the obtained total TI and mMTC slice utilities vs.
arrival rates of tactile packets. We observe that a great tactile packet
arrival rate will generally lead to small total TI and mMTC slice utilities
for IRA-ADMM and IRA. In fact, higher RU transmit power is consumed to deliver
more tactile packets, which results in a smaller TI slice utility. IRA-ADMM
obtains greater total TI and mMTC slice utilities than IRA because IRA-ADMM is
efficient than IRA in terms of bandwidth allocation. For IRA-ADMM, it
allocates bandwidth to slices in the sense of global consensus, while IRA
allocates bandwidth in a local sense.

##  IV Open issues

###  IV-A Edge intelligence

Owing to the finite speed of light, the stringent responsive requirement for
TI between remote human-robot pairs is difficult to be accommodated. To
address this issue, tactile edge intelligence schemes that promise real-time
insights by collecting and analyzing data at the tactile edge itself must be
investigated. For example, edge intelligence enables the forecast of tactile
feedback. Hence, tactile feedback can be received in real time, and then the
low-latency requirement of remote human-robot communications can be satisfied
without being limited by the human-robot distance. Besides, by predictively
caching tactile actuation, TI can also overcome the fundamental distance
limitation of tactile applications. However, to experience interactions with
physical or virtual remote objects as intuitive and natural, accurate edge
forecast results must be provided, which is challenging.

###  IV-B Fast slice reconfiguration

The dynamic characteristics (e.g., arrival and departure) of tactile traffic
will hit the resource availability in TI slices. Unfortunately, the exact
dynamic characteristics are hard to be predicted. To alleviate the effect of
prediction errors on the QoS of TI slices, slice reconfiguration is desired.
Nevertheless, the slice reconfiguration indicates some higher-layer operations
such as flow re-routing and VNF instance migration, which is time-consuming
and needs additional resources. Besides, it may cause service interruption,
which in turn degrades TI slices’ QoS. To tackle this issue, fast slice
reconfiguration with a small reconfiguration cost is worthy of study.

###  IV-C Private TI network

Safety and privacy are also regarded as crucial requirements for TI. However,
the application multiplexing on a common network may result in tactile
information exposure. To avoid information exposure, the common network is
required to have enhanced security capabilities as it generally takes a
monolithic view towards the security architecture. In this regard, a private
TI network with complete control over every aspect of the network can be
designed over the common network to transport tactile information securely and
privately. Specifically, the private TI network will be physically separate
from the common network but will share partial infrastructure of the common
network. A roaming agreement can be established between the private TI network
and the common network such that tactile information can traverse through the
common network. Certainly, the research on the private TI network is in its
infancy, and many issues such as the network architecture, the design of
flexible PHY and protocols, and the flexible spectrum sharing of the private
TI network are waiting to be investigated such that optimized services can be
provided.

##  V Conclusion

An overview of the fundamental features of CoMP-enabled RAN slicing for TI was
presented, followed by the TI architecture and the functional split of a CoMP-
enabled RAN. By integrating NFV and CU-DU network structure, the
implementation prototype of CoMP-enabled RAN slicing for TI-involved
application multiplexing was elaborated. Based on the prototype, a use case of
application multiplexing was also studied. Nevertheless, the research on CoMP-
enabled RAN slicing for TI is at its nascent stage. Some key technological
breakthroughs regarding security and low latency should be made before it can
contribute to “Tactile Internet” standards.

## References

  * [1] M. Simsek, A. Aijaz, M. Dohler, J. Sachs, and G. Fettweis, “5G-enabled Tactile Internet,” _IEEE Journal on Selected Areas in Communications_ , vol. 34, no. 3, pp. 460–473, 2016.
  * [2] A. Aijaz, M. Dohler, A. H. Aghvami, V. Friderikos, and M. Frodigh, “Realizing the Tactile Internet: Haptic communications over next generation 5G cellular networks,” _IEEE Wireless Communications_ , vol. 24, no. 2, pp. 82–89, 2016\.
  * [3] M. Series, “IMT Vision–Framework and overall objectives of the future development of IMT for 2020 and beyond,” Recommendation ITU, Tech. Rep. M.2083, Sep. 2015, https://www.itu.int/md/R12-WP5D.AR-C-0486.
  * [4] P. Rost, C. Mannweiler, D. S. Michalopoulos, C. Sartori, V. Sciancalepore, N. Sastry, O. Holland, S. Tayade, B. Han, D. Bega _et al._ , “Network slicing to enable scalability and flexibility in 5G mobile networks,” _IEEE Communications magazine_ , vol. 55, no. 5, pp. 72–79, 2017.
  * [5] 3GPP, “Coordinated multi-point operation for LTE physical layer aspects (release 11),” 3GPP, Tech. Rep. 36.819, Mar. 2011.
  * [6] N. Promwongsa, M. Abu-Lebdeh, S. Kianpisheh, F. Belqasmi, R. H. Glitho, H. Elbiaze, N. Crespi, and O. Alfandi, “Ensuring reliability and low cost when using a parallel VNF processing approach to embed delay-constrained slices,” _IEEE Transactions on Network and Service Management_ , vol. 17, no. 4, pp. 2226–2241, 2020.
  * [7] A. Ksentini, P. A. Frangoudis, P. Amogh, and N. Nikaein, “Providing low latency guarantees for slicing-ready 5G systems via two-level MAC scheduling,” _IEEE Network_ , vol. 32, no. 6, pp. 116–123, 2018.
  * [8] A. S. Shafigh, S. Glisic, and B. Lorenzo, “Dynamic network slicing for flexible radio access in Tactile Internet,” in _GLOBECOM 2017-2017 IEEE Global Communications Conference_. IEEE, 2017, pp. 1–7.
  * [9] 3GPP, “Study on new radio access technology: Radio access architecture and interfaces (release 14),” 3GPP, Tech. Rep. 38.801 V2.0.0, Apr. 2017.
  * [10] P. Yang, X. Xi, Y. Fu, T. Q. Quek, X. Cao, and D. Wu, “Multicast eMBB and bursty URLLC service multiplexing in a CoMP-enabled RAN,” _IEEE Transactions on Wireless Communications_ , 2021, in Press. DOI: 10.1109/TWC.2020.3047263.
  * [11] A. Ericsson, “Common public radio interface (CPRI),” www.cpri.info/downloads/CPRI_v_7_0_2015-10-09.pdf, 2015.
  * [12] A. Ksentini and N. Nikaein, “Toward enforcing network slicing on RAN: Flexibility and resources abstraction,” _IEEE Communications Magazine_ , vol. 55, no. 6, pp. 102–108, 2017.
  * [13] P. Yang, X. Xi, T. Q. Quek, J. Chen, X. Cao, and D. Wu, “RAN slicing for massive IoT and bursty URLLC service multiplexing: Analysis and optimization,” _IEEE Internet of Things Journal_ , 2021, in Press. DOI: 10.1109/JIOT.2021.3068518.
  * [14] C. Li, C.-P. Li, K. Hosseini, S. B. Lee, J. Jiang, W. Chen, G. Horn, T. Ji, J. E. Smee, and J. Li, “5G-based systems design for Tactile Internet,” _Proceedings of the IEEE_ , vol. 107, no. 2, pp. 307–324, 2018.

Peng Yang [M]<EMAIL_ADDRESS>is with the Singapore University of
Technology and Design (SUTD). His research interests include network slicing,
Tactile Internet, and UAV networks. Xing Xi<EMAIL_ADDRESS>is currently
pursuing the Ph.D. degree with Beihang university, Beijing, China. His
research interests include intelligent transportation systems and UAV
networks. Tony Q. S. Quek [F]<EMAIL_ADDRESS>is the Cheng Tsang Man
Chair Professor with Singapore University of Technology and Design (SUTD). He
also serves as the Head of ISTD Pillar, Sector Lead of the SUTD AI Program,
and the Deputy Director of the SUTD-ZJU IDEA. His current research topics
include wireless communications and networking, network intelligence,
internet-of-things, URLLC, and big data processing. Hyundong Shin [F]
<EMAIL_ADDRESS>is a Professor at the Department of Electronics and
Information Convergence Engineering and the Department of Electronic
Engineering, Kyung Hee University. His research interests include quantum
information science, wireless communications, and machine intelligence.
---


# QED effects on the nuclear magnetic shielding of 3He

Dominik Wehrli  Laboratorium für Physikalische Chemie, ETH-Zürich, 8093
Zürich, Switzerland  Anna Spyszkiewicz-Kaczmarek  Faculty of Chemistry, Adam
Mickiewicz University, Uniwersytetu Poznańskiego 8, 61-614 Poznań, Poland
Mariusz Puchalski  Faculty of Chemistry, Adam Mickiewicz University,
Uniwersytetu Poznańskiego 8, 61-614 Poznań, Poland  Krzysztof Pachucki
Faculty of Physics, University of Warsaw, Pasteura 5, 02-093 Warsaw, Poland

###### Abstract

The leading quantum electrodynamic corrections to the nuclear magnetic
shielding in one- and two-electron atomic systems are obtained in a complete
form, and the shielding constants of 1H, 3He+, and 3He are calculated to be
$17.735\,436(3)\cdot 10^{-6}$, $35.507\,434(9)\cdot 10^{-6}$, and
$59.967\,029(23)\cdot 10^{-6}$, respectively. These results are orders of
magnitude more accurate than previous ones, and, with the ongoing measurement
of the nuclear magnetic moment of 3He+ and planned 3He2+, they open the window
for high-precision absolute magnetometry using 3He NMR probes. The presented
theoretical approach is applicable to all other light atomic and molecular
systems, which facilitates the improved determination of magnetic moments of
any light nuclei.

††preprint: Version 1.0

Introduction – The chemical inertness of 3He atomic gas and the absence of an
electric quadrupole moment yield an exceptional isolation of nuclear spins
that form a macroscopic quantum state achieving long coherence times. Because
it can be hyper-polarized using laser optical pumping techniques, resulting in
a very stable and sensitive probe to the magnetic field, the nuclear magnetic
resonance (NMR) magnetometry founded on 3He is of great interest for different
areas of physics and applied sciences [1]. In particular, 3He gas cells have
recently been used in several fundamental physics experiments for their
capability of determining absolute field values, in the calibration of
magnetic sensors, e.g., in muon $g-2$ measurements [2], in searching for a
permanent electric dipole moment [3], and in the development of cryogenic NMR
techniques for new experiments on the electron and positron magnetic moments
[4, 5]. For this reason, hyper-polarized 3He NMR probes have been proposed as
a new standard for absolute magnetometry [6, 7, 5], and consequently, a high-
precision value of the helion nuclear magnetic moment $\mu_{h}$ would be
indispensable.

For the determination of $\mu_{h}$, the NMR measurements of 3He performed with
respect to H2 and H2O [8, 9, 10, 11, 12] have been used to obtain the current
CODATA value [13] with a relative accuracy of $1.2\cdot 10^{-8}$. This
approach, based on a comparison to the proton magnetic moment $\mu_{p}$, is
limited by the little-known magnetic shielding effects caused by the
surrounding particles (electrons, nuclei). To obtain accurate shielding
factors, we can consider using theoretical methods. Nonetheless, this
possibility only applies to small atomic and molecular systems, for which we
can use an approach based on nonrelativistic implementation of quantum
electrodynamics (NRQED). Such calculations in H2 with an accuracy similar to
the magnetic shielding of 3He performed in this work can lead to $\mu_{h}$
with relative accuracy as high as of that of $\mu_{p}$, i.e., 2.9 $\cdot
10^{-10}$ [14].

In another approach, the 3He magnetic moment is related to that of the
electron. This is achieved by measuring the ratio of the magnetic moment of
the $2^{3}S$ metastable state of 4He to the ground state of 3He [15], and,
independently, by a very recent measurement of the Breit-Rabi splitting of
3He+ in a Penning trap [16, 17]. This latter experiment may provide an
improved result for $\mu_{h}$ in comparison to the current CODATA value.

Finally, an effort is underway [16, 14] to directly measure $\mu_{h}$ in a
cryogenic Penning trap using techniques similar to those applied for the
proton [18, 19] and the antiproton [20]. The high-precision results of the
magnetic shielding constant obtained in this work provide the shielded
magnetic moment of 3He+ and 3He at the accuracy of the measured value of
$\mu_{h}$, which is a prerequisite for the realization of 3He NMR probes as a
new standard for absolute magnetometry.

Magnetic shielding – Let us recall the definition of the nuclear diamagnetic
shielding of a single atom [21, 22]. The coupling of the nuclear magnetic
moment $\vec{\mu}$ with the homogeneous magnetic field $\vec{B}$ modified by
the presence of atomic electrons can be parametrized in terms of the nuclear
magnetic shielding constant $\sigma$ according to

$\delta H=-\vec{\mu}\cdot\vec{B}\,(1-\sigma).$ (1)

In order to calculate accurately the parameter $\sigma$ we employ the so-
called nonrelativistic QED [23], and assume the expansion of $\sigma$ as a
double power series in the fine structure constant $\alpha$ and the electron-
nucleus mass ratio $m/m_{N}$,

$\displaystyle\sigma$
$\displaystyle=\sigma^{(2)}+\sigma^{(4)}+\sigma^{(5)}+\sigma^{(6)}+\sigma^{(2,1)}+\sigma^{(2,2)}+\sigma^{(4,1)}+\ldots$
(2)

In this equation, $\sigma^{(n)}\propto\alpha^{n}$ coefficients are the
nonrelativistic shielding, the relativistic, the leading QED, and the higher-
order QED corrections in the infinite nuclear mass approximation,
respectively. The terms $\sigma^{(n,k)}\propto\alpha^{n}\,(m/m_{N})^{k}$ are
corrections due to the finite nuclear mass. The main advantage of this
expansion is the possibility to derive exact formulas for expansion
coefficients in terms of some matrix element with the nonrelativistic wave
function, which can be accurately calculated.

In the case of one- and two-electron atomic systems all lower-order
contributions in Eq. (2) are known in the literature. Namely, the formula for
the leading term $\sigma^{(2)}$ valid for atoms

$\displaystyle\sigma^{(2)}=\frac{\alpha}{3\,m}\,\Bigl{\langle}\sum_{a}\frac{1}{r_{a}}\Bigr{\rangle}$
(3)

was first introduced by Lamb [21] and later generalized for molecules by
Ramsey [22]. The correction $\sigma^{(4)}$ can be found in Refs. [24, 25] and
references therein, while $\sigma^{(2,1)}$ and $\sigma^{(2,2)}$ were derived
in Ref. [26]. In contrast, QED corrections until now have only been partially
investigated and therefore are presently the bottleneck of theoretical
predictions. There have been attempts [27] to include them within the
formalism based on the Dirac-Coulomb (DC) Hamiltonian [28], but there is
currently no adequate formulation of the QED theory for many electron systems.
While, for hydrogenic ions, Yerokhin et al. [29, 30] performed a
nonperturbative numerical evaluation of one-loop QED contributions and
observed a slow numerical convergence for the small nuclear charge $Z$.
Therefore, these results were supplemented by direct NRQED evaluation of the
leading QED correction $\sim\alpha^{5}$. However, some effects due to the
magnetic moment anomaly were omitted there, which was the reason for small
discrepancies with the nonperturbative results for the medium-Z hydrogen-like
ions. For helium, only the leading QED logarithmic correction
$\sim\alpha^{5}\,\ln\alpha$ was obtained by Rudziński et al. in Ref. [25]. In
this Letter, we present the complete formulas for $\sigma^{(5)}$, valid for
one- and two-electron atomic systems, which can be easily generalized to any
light few-electron system. They are supplemented by the numerical calculations
performed for 1H, 3He+ and 3He, using explicitly correlated exponential
functions for the latter case, what ensures high precision of numerical
results. We will use $\hbar=c=\varepsilon_{0}=1$ throughout the paper and the
CODATA 2018 values of physical constants [13].

NRQED approach – In the NRQED formalism we can include QED effects, coming
from large photon momenta, in the generalized Breit-Pauli Hamiltonian (see
e.g. [31]). For the case of a two-electron atomic system coupled to a magnetic
field $\vec{B}$, this effective Hamiltonian is given by

$H_{\mathrm{BP}}=\sum_{a=1}^{2}H_{a}+H_{12},$ (4) $\displaystyle H_{a}$
$\displaystyle=\frac{\vec{\pi}_{a}^{2}}{2\,m}-\frac{\vec{\pi}_{a}^{4}}{8\,m^{3}}-\frac{Z\,\alpha}{r_{a}}-\frac{e\,(1+\kappa)}{2\,m}\,\vec{\sigma}_{a}\cdot\vec{B}_{a}+\frac{e}{8\,m^{3}}\,\\{\vec{\pi}_{a}^{2},\vec{\sigma}_{a}\cdot\vec{B}_{a}\\}-\frac{e^{2}}{2}\bigg{(}\frac{1}{4\,m^{3}}+\alpha_{M}\bigg{)}\vec{B}_{a}^{2}+\frac{e\,\kappa}{8\,m^{3}}\,\\{\vec{\pi}_{a}\cdot\vec{B}_{a},\vec{\sigma}_{a}\cdot\vec{\pi}_{a}\\}$
$\displaystyle+\frac{Z\,\alpha\,(1+2\,\kappa)}{4\,m^{2}}\,\vec{\sigma}_{a}\cdot\frac{\vec{r}_{a}\times\vec{\pi}_{a}}{r_{a}^{3}}+\frac{e}{6\,m}\,\bigg{(}r^{2}_{E}+r^{2}_{\rm
vp}-\frac{3\,\kappa}{4\,m^{2}}\bigg{)}(\nabla_{a}\times\vec{B}_{a})\cdot\vec{\pi}_{a}+\frac{2\,\pi\,Z\,\alpha}{3}\bigg{(}\frac{3}{4\,m^{2}}+r^{2}_{E}+r^{2}_{\rm
vp}\bigg{)}\,\delta(\vec{r_{a}}),$ (5) $\displaystyle H_{12}$
$\displaystyle=\frac{\alpha}{r}-\frac{4\,\pi\,\alpha}{3}\left(\frac{3}{4\,m^{2}}+r_{E}^{2}+\frac{1}{2}\,r^{2}_{\rm
vp}+\frac{\,(1+\kappa)^{2}}{2\,m^{2}}\,\vec{\sigma}_{1}\cdot\vec{\sigma}_{2}\right)\delta(\vec{r})+\frac{\alpha}{4\,m^{2}\,r^{3}}\Big{[}2\,(1+\kappa)\big{(}\vec{\sigma}_{1}\cdot\vec{r}\times\vec{\pi}_{2}-\vec{\sigma}_{2}\cdot\vec{r}\times\vec{\pi}_{1}\big{)}$
$\displaystyle\
+(1+2\,\kappa)\big{(}\vec{\sigma}_{2}\cdot\vec{r}\times\vec{\pi}_{2}-\vec{\sigma}_{1}\cdot\vec{r}\times\vec{\pi}_{1}\big{)}\Big{]}-\frac{\alpha}{2\,m^{2}}\,\pi_{1}^{i}\left(\frac{\delta^{ij}}{r}+\frac{r^{i}\,r^{j}}{r^{3}}\right)\pi_{2}^{j}-\frac{3\,\alpha\,(1+\kappa)^{2}}{4\,m^{2}}\,\frac{\sigma_{1}^{i}\,\sigma_{2}^{j}}{r^{3}}\biggl{(}\frac{r_{1}^{i}\,r_{2}^{j}}{r^{2}}-\frac{\delta^{ij}}{3}\biggr{)}$
(6)

where $\vec{r}=\vec{r}_{1}-\vec{r}_{2}$, $\vec{\pi}=\vec{p}-e\,\vec{A}$ is the
generalized momentum in the external field, $\kappa=\alpha/(2\,\pi)$ is the
magnetic moment anomaly, and

$\displaystyle r_{E}^{2}$
$\displaystyle=\frac{3\,\kappa}{2\,m^{2}}+6\,F^{\prime}_{1}(0)=\frac{2\,\alpha}{\pi\,m^{2}}\left(\ln\frac{m}{2\,\epsilon}+\frac{5}{6}\right),$
(7) $\displaystyle r^{2}_{\rm vp}$
$\displaystyle=-\frac{2\,\alpha}{5\,\pi\,m^{2}},$ (8)
$\displaystyle\alpha_{M}$
$\displaystyle=\frac{4\,\alpha}{3\,\pi\,m^{3}}\left(-\ln\frac{m}{2\,\epsilon}+\frac{13}{24}\right).$
(9)

The parameters $r_{E}$ and $\alpha_{M}$ are interpreted as the charge radius
and the magnetic polarizability of an electron, respectively. They depend on
$\epsilon$, i.e. the photon momentum cut-off being used as a regulator [32],
and this dependence cancels out in the complete expression for any physical
quantity. The formula for $r^{2}_{E}$ is derived from the known radiative
correction to electromagnetic form factors $F_{1}$ and $F_{2}$ [32], and
$r^{2}_{\rm vp}$ incorporates the corrections due to the vacuum polarization.
The formula for $\alpha_{M}$ has for the first time been presented in Ref.
[30] and can be obtained in a similar way as for the electric polarizability,
denoted by $\chi$, in Ref. [33]. Without QED ($r^{2}_{E}=r^{2}_{\rm
vp}=\alpha_{M}=\kappa=0$) $H_{\mathrm{BP}}$ is the standard Breit-Pauli
Hamiltonian [34].

For the derivation of the magnetic shielding constant $\sigma$, we consider
Eq. (4) for an atomic system in a magnetic field corresponding to the sum of
the vector potential $\vec{A}_{I}$ due to the magnetic moment $\vec{\mu}$ of
the nucleus , and $\vec{A}_{E}$ due to the homogeneous external magnetic field
$\vec{B}$, namely

$\vec{A}_{E}=\frac{1}{2}\,\vec{B}\times\vec{r},\quad{\rm
and}\quad\vec{A}_{I}=\frac{1}{4\,\pi}\,\vec{\mu}\times\frac{\vec{r}}{r^{3}}.$
(10)

Following Ramsey’s theory of the magnetic shielding [22, 35], we split the
Hamiltonian $H_{\mathrm{BP}}$ as

$\displaystyle H_{\mathrm{BP}}$ $\displaystyle=H_{0}+\delta
H_{\vec{A}_{E}=\vec{A}_{I}=0}+\delta H_{\vec{A}_{I},\vec{A}_{E}=0}+\delta
H_{\vec{A}_{E},\vec{A}_{I}=0}$ $\displaystyle+\delta
H_{\vec{A}_{E},\vec{A}_{I}}+O(\vec{A}^{2}_{I,E}),$ (11)

where $\delta H$ is treated as a perturbation to the nonrelativistic
Hamiltonian $H_{0}$, $\delta H_{\vec{A}_{E}=\vec{A}_{I}=0}$ is independent of
the magnetic fields, $\delta H_{\vec{A}_{I},\vec{A}_{E}=0}$ is linear in
$\vec{A}_{I}$, $\delta H_{\vec{A}_{E},\vec{A}_{I}=0}$ is linear in
$\vec{A}_{E}$, and $\delta H_{\vec{A}_{E},\vec{A}_{I}}$ is bilinear in both
fields. Because we are only interested in energy corrections that are
proportional to $\vec{\mu}\cdot\vec{B}$, we write

$\displaystyle\delta E$ $\displaystyle=\langle\delta
H_{\vec{A}_{E},\vec{A}_{I}}\rangle+2\left\langle\delta
H_{\vec{A}_{E},\vec{A}_{I}}\,\frac{1}{(E_{0}-H_{0})^{\prime}}\,\delta
H_{\vec{A}_{E}=\vec{A}_{I}=0}\right\rangle$ $\displaystyle+2\left\langle\delta
H_{\vec{A}_{I},\vec{A}_{E}=0}\,\frac{1}{(E_{0}-H_{0})^{\prime}}\,\delta
H_{\vec{A}_{E},\vec{A}_{I}=0}\right\rangle+\ldots\,,$ (12)

where $1/(E_{0}-H_{0})^{\prime}$ is the reduced Green’s function, and the
ellipses denote terms that are not proportional to $\vec{\mu}\cdot\vec{B}$ and
will be discarded. The expectation values are taken with respect to the ground
electronic state of $H_{0}$, which is an $S$ state in the case of hydrogen-
and helium-like systems. The spherical symmetry then implies the relation
$\mu^{i}B^{j}=\delta^{ij}/3\,\vec{\mu}\cdot\vec{B}$ allowing for a simple
factoring of $\vec{\mu}\cdot\vec{B}$ from many terms appearing in Eq. (12),
and $\sigma$ is obtained through the relation $\delta
E=\vec{\mu}\cdot\vec{B}\,\sigma$.

Leading QED correcion – We derive the correction $\sigma^{(5)}$ for a helium-
like system, bearing in mind the corresponding derivation of the Lamb shift
[36]. Namely, $\sigma^{(5)}$ is given by the sum of the high $\sigma_{B}$ and
the low $\sigma_{A}$ energy parts

$\sigma^{(5)}=\sigma_{B0}+\sigma_{B1}+\sigma_{B2}+\sigma_{B3}+\sigma_{A1}+\sigma_{A2}.$
(13)

The high energy part is obtained as follows. $\sigma_{B0}$ corresponds to the
helium Lamb shift [36] with the wave function corrected by the effect of
magnetic fields in Eq. (3)

$\displaystyle\sigma_{B0}$
$\displaystyle=\frac{2\,\alpha^{3}}{3\,m^{3}}\,\biggl{\langle}{\mkern-3.0mu}\biggl{(}\frac{1}{r_{1}}+\frac{1}{r_{2}}\biggr{)}\frac{1}{(E_{0}-H_{0})^{\prime}}\biggl{\\{}{\mkern-3.0mu}\left(\frac{19}{30}+\ln(\alpha^{-2})\right)$
$\displaystyle\times\frac{4\,Z}{3}\,(\delta(\vec{r_{1}})+\delta(\vec{r_{2}}))+\left(\frac{164}{15}+\frac{14}{3}\,\ln\alpha\right)\delta(\vec{r})$
$\displaystyle-\frac{7\,\alpha^{3}\,m^{3}}{6\,\pi}\,P\left(\frac{1}{(m\,\alpha\,r)^{3}}\right){\mkern-3.0mu}\biggr{\\}}{\mkern-3.0mu}\biggr{\rangle}.$
(14)

where $P(1/r^{3})$ is the Araki-Sucher term [36]. The next parts are directly
obtained from the NRQED Hamiltonian in Eq. (4), and are the following

$\displaystyle\sigma_{B1}$
$\displaystyle=\frac{\alpha^{2}}{m^{3}}\left(\frac{20}{9}\,\ln(\alpha^{-2})-\frac{1361}{540}\right)\left\langle\delta(\vec{r}_{1})+\delta(\vec{r}_{2})\right\rangle,$
(15) $\displaystyle\sigma_{B2}$
$\displaystyle=-\frac{\alpha^{2}}{3\,m^{3}}\bigg{\langle}{\mkern-5.0mu}(\delta(\vec{r}_{1})-\delta(\vec{r}_{2}))\frac{1}{(E_{0}-H_{0})}$
(16)
$\displaystyle\times\bigg{(}\frac{4\,\vec{p}_{1}^{\;2}}{3\,m}-\frac{4\,\vec{p}_{2}^{\;2}}{3\,m}-\frac{Z\,\alpha}{r_{1}}+\frac{Z\,\alpha}{r_{2}}-\frac{\alpha}{3\,r^{3}}\,\vec{r}\cdot(\vec{r}_{1}+\vec{r}_{2})\bigg{)}{\mkern-5.0mu}\bigg{\rangle},$
$\displaystyle\sigma_{B3}$
$\displaystyle=-\frac{3\,\alpha^{2}}{16\,\pi\,m^{3}}\bigg{\langle}{\mkern-5.0mu}\bigg{(}\frac{r_{1}^{i}\,r_{1}^{j}}{r_{1}^{5}}-\frac{r_{2}^{i}\,r_{2}^{j}}{r_{2}^{5}}\bigg{)}^{(2)}{\mkern-20.0mu}\frac{1}{(E_{0}-H_{0})}\bigg{(}Z\,\alpha\,\frac{r_{1}^{i}\,r_{1}^{j}}{r_{1}^{3}}$
$\displaystyle-Z\,\alpha\,\frac{r_{2}^{i}\,r_{2}^{j}}{r_{2}^{3}}+\alpha\,\frac{r^{i}\,(r_{1}^{j}+r_{2}^{j})}{3\,r^{3}}+\frac{2}{3\,m}\,(p_{1}^{i}\,p_{1}^{j}-p_{2}^{i}\,p_{2}^{j})\bigg{)}{\mkern-5.0mu}\bigg{\rangle}.$
(17)

where the second rank tensor is defined by
$(p^{i}q^{j})^{(2)}\equiv(p^{i}q^{j}+p^{j}q^{i})/2-\delta^{ij}/3\,\vec{p}\cdot\vec{q}$.

Analogously to the calculation of the Lamb shift [36], we define the low-
energy energy contribution as

$\delta
E_{A}=-\frac{2\,\alpha}{3\,\pi}\left\langle{\mkern-5.0mu}(\vec{\pi}_{1}+\vec{\pi}_{2})\,(H-E)\ln\frac{2\,(H-E)}{m\,\alpha^{2}}\,(\vec{\pi}_{1}+\vec{\pi}_{2}){\mkern-5.0mu}\right\rangle_{{\mkern-6.0mu}B}{\mkern-6.0mu},$
(18)

where $\langle\ldots\rangle_{B}$ denotes the expectation value with respect to
the ground state with energy $E$ of the Hamiltonian in the presence of the
magnetic field

$\displaystyle H$
$\displaystyle=\frac{\vec{\pi}^{2}_{1}}{2\,m}+\frac{\vec{\pi}^{2}_{2}}{2\,m}-\frac{Z\,\alpha}{r_{1}}-\frac{Z\,\alpha}{r_{2}}+\frac{\alpha}{r}$
$\displaystyle=H_{0}+\frac{1}{3\,m}\,\vec{\mu}\cdot\vec{B}\,U-\frac{e}{2\,m}\,\vec{L}\cdot\vec{B}-\frac{e}{4\,\pi\,m}\vec{\mu}\cdot\vec{U},$
(19)

where $L^{i}=L_{1}^{i}+L_{2}^{i}$, $U=\alpha/r_{1}+\alpha/r_{2}$, and
$U^{i}=L_{1}^{i}/r_{1}^{3}+L_{2}^{i}/r_{2}^{3}$. Equation (18) is only a
formal expression for $\delta E_{A}$, and it needs to be expanded in the
magnetic field. For this we rewrite $\delta E_{A}$ in the form

$\delta
E_{A}=-\frac{2\,\alpha}{3\,\pi}\left\langle{\mkern-5.0mu}(\vec{r}_{1}+\vec{r}_{2})\,(H-E)^{3}\ln\frac{2\,(H-E)}{m\,\alpha^{2}}\,(\vec{r}_{1}+\vec{r}_{2}){\mkern-5.0mu}\right\rangle_{{\mkern-6.0mu}B}{\mkern-6.0mu},$
(20)

and divide $\delta E_{A}$ into two parts coming from different perturbations
to $H_{0}$ in Eq. (19). The part due to $1/(3\,m)\,\vec{\mu}\cdot\vec{B}\,U$
leads to the first correction

$\displaystyle\sigma_{A1}$
$\displaystyle=-\frac{2\,\alpha}{9\,\pi}\,\delta_{U}\bigg{\langle}{\mkern-5.0mu}(\vec{r}_{1}+\vec{r}_{2})\,(H_{0}-E_{0})^{3}\ln\frac{2\,(H_{0}-E_{0})}{m\,\alpha^{2}}$
$\displaystyle\times(\vec{r}_{1}+\vec{r}_{2}){\mkern-5.0mu}\bigg{\rangle}=\frac{2\,\alpha}{9\,\pi\,m^{2}}\,{\cal
D}_{A1}\,\ln k_{A1},$ (21) $\displaystyle{\cal D}_{A1}$ $\displaystyle\equiv\
-4\,\pi\,Z\bigg{\langle}U\,\frac{1}{(E_{0}-H_{0})^{\prime}}\,(\delta(\vec{r_{1}})+\delta(\vec{r_{2}}))\bigg{\rangle}$
$\displaystyle+2\,\pi\,\alpha\langle\delta(\vec{r_{1}})+\delta(\vec{r_{2}})\rangle$
(22)

The second correction is due to perturbation from the two other terms in Eq.
(19),

$\displaystyle\sigma_{A2}$
$\displaystyle=-\frac{\alpha^{2}}{9\,\pi\,m^{2}}\delta_{L^{i},U^{i}}\bigg{\langle}{\mkern-5.0mu}(\vec{r}_{1}+\vec{r}_{2})(H_{0}-E_{0})^{3}\ln\frac{2\,(H_{0}-E_{0})}{m\,\alpha^{2}}$
$\displaystyle\times(\vec{r}_{1}+\vec{r}_{2}){\mkern-5.0mu}\bigg{\rangle}=-\frac{\alpha^{2}}{9\,\pi\,m^{2}}\,(1+3\,\ln
k_{A2})\,{\cal D}_{A2},$ (23) $\displaystyle{\cal D}_{A2}$
$\displaystyle\equiv
8\,\pi\langle\delta(\vec{r_{1}})+\delta(\vec{r_{2}})\rangle$ (24)

Both $\sigma_{A1}$ and $\sigma_{A2}$ are expressed with the Bethe-type
logarithms $\ln k_{A1}$ and $\ln k_{A2}$, respectively. Their numerical
calculation is not straightforward because it involves the logarithm of the
Hamiltonian (see e.g. [37]). Like the standard Bethe logarithm $\ln k_{0}$,
all these elements have the striking property that they are only slightly
dependent on the number of electrons, as can be seen from the numerical
results in the Table 1.

Hydrogenic formula – The expression $\sigma^{(5)}$ for a one-electron atomic
system is obtained by direct reduction of one of the electrons in Eq. (13) and
related expressions. As a result, the total magnetic shielding for hydrogen-
like ions, including contributions up to order $\alpha^{5}$, has a compact
structure, namely

$\displaystyle\sigma$
$\displaystyle=\frac{1}{3}\,\alpha\,(Z\,\alpha)+\frac{97}{108}\,\alpha\,(Z\,\alpha)^{3}+\frac{8}{9\,\pi}\,\alpha^{2}\,(Z\,\alpha)^{3}\,\biggl{[}\ln\left[(Z\,\alpha)^{-2}\right]$
$\displaystyle+2\ln k_{0}(H)-3\ln
k_{3}-\frac{221}{64}+\frac{3}{5}\biggr{]}+\frac{Z\,\alpha^{2}}{3}\,\biggl{(}\frac{1-g_{N}}{g_{N}}-1\biggr{)}\frac{m}{m_{N}}$
$\displaystyle+\frac{Z\,\alpha^{2}}{3}\left(1+Z\,\frac{2+3\,g_{N}}{2\,g_{N}}-2\,\frac{1-g_{N}}{g_{N}}\right)\frac{m^{2}}{m_{N}^{2}}.$
(25)

where the nuclear $g$-factor is defined as

$g_{N}=\frac{m_{N}}{Z\,m_{p}}\,\frac{\mu}{\mu_{N}}\,\frac{1}{I}.$ (26)

In Eq. (26), $m_{p}$ is the proton mass, $\mu_{N}$ is the nuclear magneton,
and $\mu$ and $I$ are the magnetic moment and the spin of the considered
nucleus, respectively. Moreover, the $\ln k_{A1}$ and $\ln k_{A2}$ introduced
in Eqs. (21) and (23) can be represented for hydrogenic systems in terms of
$\ln k_{0}({\rm H})$ and $\ln k_{3}$ [38] (as shown in Table 1).

Comparing Eq. (25) to the result published in Refs. [29, 30] it differs in the
constant term $-421/96$ instead of $-221/64$, due to omitted contributions
coming from magnetic moment anomaly. We observe that the new analytic result
for hydrogenic systems is in much better agreement with the numerical values
from Refs. [29, 30] which were calculated to all orders in $Z\,\alpha$ but
exhibited large numerical uncertainties for $Z<10$.

Numerical calculations and results – To evaluate the magnetic shielding
constant for 3He including the QED correction $\sigma^{(5)}$ we represent the
wave function and other auxiliary functions used in the calculation of the
second-order terms in the basis set of explicitly correlated exponential
functions. For example, the ground state wave function is given by

$\psi(\vec{r}_{1},\vec{r}_{2})=(1+{\cal
P}_{12})\sum_{i=1}^{N}c_{i}\,e^{-\alpha_{i}\,r_{1}-\beta_{i}\,r_{2}-\gamma_{i}\,r_{12}}\,,$
(27)

where ${\cal P}_{12}$ exchanges $\vec{r}_{1}$ with $\vec{r}_{2}$, and all
nonlinear parameters $\alpha_{i}$, $\beta_{i}$, and $\gamma_{i}$ were
determined variationally. In order to control the numerical uncertainty, we
performed the calculations with several basis sets successively increasing
their size by a factor of two; e.g., for the ground-state wave function we
used $N=128$, $256$, and $512$ basis functions. The result obtained in the
largest basis set is accurate to 17 significant digits compared to the
benchmark result obtained by Korobov at al. [39]. The method of full
optimization used in this work has some advantages in numerical calculations.
First of all, the numerical uncertainties with these compact wave functions
are negligible in comparison to uncertainties due to the omitted higher-order
corrections. Moreover, the calculations of all the second-order matrix
elements have been performed in quadrupole precision arithmetics, with the
full optimization of the intermediate basis of the comparable same size as the
ground-state wave function. Similarly, the optimization of the intermediate
basis in the calculation of Bethe-type logarithms using integral
representations was essential in obtaining high accuracy of numerical results.

Table 1: Numerical values of He matrix elements. The results are in atomic units and using the notation $\sigma^{(n)}=\alpha^{n}\tilde{\sigma}^{(n)}$.  operator | value
---|---
$E_{0}$ | $-2.903\,724\,377\,034\,119\,59(1)$
$\left\langle U\right\rangle$ | $3.376\,633\,601\,434\,081(4)$
${\cal D}_{A1}$ | $98.798\,613\,9(3)$
${\cal D}_{A2}$ | $91.002\,103\,1(3)$
$\ln k_{0}$ | $4.370\,160\,22(2)$
$\ln k_{0}({\rm H})+2\,\ln Z$ | $4.370\,422\,917$
$\ln k_{A1}$ | $4.829\,409(3)$
$\ln k_{0}({\rm H})+2\,\ln Z+1/2$ | $4.870\,422\,917$
$\ln k_{A2}$ | $4.638\,660(15)$
$\ln k_{3}+2\,\ln Z$ | $4.659\,100\,906$
$\tilde{\sigma}_{A1}$ | $33.750\,67(2)$
$\tilde{\sigma}_{A2}$ | $-48.007\,69(14)$
$\tilde{\sigma}_{B0}$ | $70.054\,125\,1(2)$
$\tilde{\sigma}_{B1}$ | $-55.342\,119\,09(14)$
$\tilde{\sigma}_{B2}$ | $4.188\,033\,454(7)$
$\tilde{\sigma}_{B3}$ | $0.011\,67(3)$
$\tilde{\sigma}^{(5)}$ | $4.654\,69(15)$

From the analysis of convergence, we obtained the extrapolated mean values of
the operators. The different contributions to the QED correction of the
magnetic shielding in 3He are given in Table 1, while the numerical
calculation of the lower-order contributions in $\alpha$ is presented in
detail in Ref. [25]. The numerical value for the shielding Bethe logarithms
$\ln k_{A1}$ and $\ln k_{A2}$ were verified in Table 1 by the hydrogenic
counterparts, because we expect a minor dependence on the number of electrons
similarly to standard Bethe logarithms.

The final results for the shielding in 3He are summarized in Table 2. We note
that the QED correction is subject to cancellations of the different
contributions, and the value presented here for $\sigma^{(5)}=96.3\cdot
10^{-12}$ is significantly lower than the one previously reported in Ref.
[25], $\sigma^{(5)}_{\mathrm{prev}}=502\cdot 10^{-12}$, in which only the
leading logarithmic contribution was included. This difference is also
reflected in the final value of the shielding. The convergence of the
expansions in the fine structure constant $\alpha$ and the electron-nuclear
mass ratio $m/m_{N}$ is very rapid, which justifies our approach based on the
NRQED theory. The numerical uncertainties are completely negligible, and only
the unknown higher-order terms in $\alpha$ and $m/m_{N}$ contribute to the
uncertainty. The finite nuclear size effects $E_{\rm NS}$ are significant only
for heavy elements; here we expect them to be as important as for the binding
energy, which for the ground state of He amounts to [40] $E_{\rm
NS}/E_{0}\approx 5\cdot 10^{-9}$, which is negligible with respect to the
current relative uncertainty of $\sigma$.

Table 2: Contributions to the shielding constant $\sigma\cdot 10^{6}$ for 1H, 3He+, and 3He. Quantities that are preceded by “$\pm$” represent the uncertainties. Numerical values $\sigma^{(2)}$ and $\sigma^{(4)}$ for 3He are taken from Ref. [25]. $\sigma^{(2,2)}$(He) is estimated to be between $\sigma^{(2,2)}$(He+) and $2\,\sigma^{(2,2)}$(He+). The relative uncertainty of the finite nuclear mass correction $\sigma^{(4,1)}$ is estimated as $2\,m/m_{N}$ of $\sigma^{(4)}$. $\sigma^{(6)}$ is partially known from the Dirac equation [24], but we expect cancellation with the radiative correction, so we estimate the uncertainty originating from this contribution as $(Z\,\alpha)^{2}\,\sigma^{(4)}$. |  |  |
---|---|---|---
|  1H |  3He+ |  3He
$\sigma^{(2)}$ | $17.750\,451\,5$ | $35.500\,903\,0$ | $59.936\,770\,5$
$\sigma^{(4)}$ | $0.002\,546\,9$ | $0.020\,375\,1$ | $0.052\,663\,1$
$\sigma^{(5)}$ | $0.000\,018\,4$ | $0.000\,082\,0$ | $0.000\,096\,3$
$\sigma^{(2,1)}$ | $-0.017\,603\,7$ | $-0.013\,933\,4$ | $-0.022\,511\,5$
$\sigma^{(2,2)}$ | $0.000\,022\,7$ | $0.000\,007\,1$ | $0.000\,010\,7(36)$
$\sigma^{(4,1)}$ | $\pm 0.000\,002\,8$ | $\pm 0.000\,007\,4$ | $\pm 0.000\,019\,2$
$\sigma^{(6)}$ | $\pm 0.000\,000\,1$ | $\pm 0.000\,004\,3$ | $\pm 0.000\,011\,2$
$\sigma\cdot 10^{6}$ | $17.735\,436(3)$ | $35.507\,434(9)$ | $59.967\,029(23)$
Vaara, Pyykkö (2003) [41] |  | $59.93(4)$
Kudo, Fukui (2005) [42] |  | $59.8$
Antušek et al. (2007) [43] |  | $59.908\,03$
Rudziński et al. (2009) [25] |  | $59.967\,43(10)$
Seino, Hada (2010) [44] |  | $59.95$
Kupka et al. (2013) [45] |  | $59.930$
|  |  |

Our result for the magnetic shielding in 3He is in disagreement with almost
all previous calculations in Refs. [42, 43, 44, 45], which use the standard
quantum chemistry codes and do not present any uncertainties. Moreover, none
of them include the finite nuclear mass or the QED effects, but the
differences are significantly greater than the omitted corrections. The only
exception is the results of Vaara and Pyykkö [41], which, although the oldest
one, present uncertainties in agreement with our result.

Conclusions – We have presented the complete formula of the leading quantum
electrodynamics correction to the magnetic shielding for hydrogen- and helium-
like atomic systems, and performed calculations for 1H and the particularly
important cases of 3He+ and 3He, which are accurate to $3\cdot 10^{-12}$,
$9\cdot 10^{-12}$, and $23\cdot 10^{-12}$, respectively. The results permit
establishing a new standard for the absolute magnetometry based on 3He NMR
probes, by removing one of the bottlenecks in determining $\mu_{h}$ by
indirect methods.

Moreover, the NMR frequency ratio of 3He and of H2 is known with a relative
accuracy of about $10^{-9}$ [12] and $\mu_{p}$ to $2.9\cdot 10^{-10}$, so it
is enough to provide the magnetic shielding of H2 with about $10^{-5}$
relative accuracy, to obtain $\mu_{h}$ with about $10^{-9}$ accuracy for
verifying the consistency among various determinations. This would require the
calculation of only nonadiabatic and relativistic effects, which we plan to
perform this in the near future, while the QED corrections, estimated using
the atomic value in Table 2, are negligible.

The result of the magnetic shielding in 3He+ can be used for the determination
of $\mu_{h}$ from the measurement of 3He+ [16, 14], being close to completion
[17]. Alternatively, if $\mu_{h}$ is measured directly by comparison with the
cyclotron frequency [14], the result for the 3He magnetic shielding will
provide an accurate magnetic moment of atomic 3He for use in absolute
magnetometry. Finally, the NRQED approach, presented in this work, to
calculate the magnetic shielding, can be extended to nuclei of other light
atomic systems, such as Li, Be, and B, which will lead to the determination of
their magnetic moments with the accuracy of experiments, which for the case of
Be+ [46] is about $10^{-9}$, i.e. higher by a few orders of magnitude in
comparison to the presently known value of $\mu_{\rm Be}$ [47, 48].

###### Acknowledgements.

D.W. thanks F. Merkt for his unconditional support to work on this project.
This research was supported by National Science Center (Poland) Grants No.
2017/27/B/ST2/02459 and 2019/35/N/ST4/04445, as well as by a computing grant
from the Poznań Supercomputing and Networking Center and by PL-Grid
Infrastructure.

## References

  * Gentile  _et al._ [2017] T. R. Gentile, P. J. Nacher, B. Saam, and T. G. Walker, Rev. Mod. Phys. 89, 045004 (2017).
  * Abi and et al. [2021] B. Abi and et al. (Muon $g-2$ Collaboration), Phys. Rev. Lett. 126, 141801 (2021).
  * Sachdeva  _et al._ [2019] N. Sachdeva, I. Fan, E. Babcock, M. Burghoff, T. E. Chupp, S. Degenkolb, P. Fierlinger, S. Haude, E. Kraegeloh, W. Kilian, S. Knappe-Grüneberg, F. Kuchler, T. Liu, M. Marino, J. Meinel, K. Rolfs, Z. Salhi, A. Schnabel, J. T. Singh, S. Stuiber, W. A. Terrano, L. Trahms, and J. Voigt, Phys. Rev. Lett. 123, 143003 (2019).
  * Gabrielse  _et al._ [2019] G. Gabrielse, S. E. Fayer, T. G. Myers, and X. Fan, Atoms 7, 45 (2019).
  * Fan  _et al._ [2019] X. Fan, S. E. Fayer, and G. Gabrielse, Rev. Sci. Instrum. 90, 083107 (2019).
  * Nikiel  _et al._ [2014] A. Nikiel, P. Blümler, W. Heil, M. Hehn, S. Karpuk, A. Maul, E. Otten, L. M. Schreiber, and M. Terekhov, Eur. Phys. J. D 68, 330 (2014).
  * Farooq  _et al._ [2020] M. Farooq, T. Chupp, J. Grange, A. Tewsley-Booth, D. Flay, D. Kawall, N. Sachdeva, and P. Winter, Phys. Rev. Lett. 124, 223001 (2020).
  * Neronov and Seregin [2014] Y. I. Neronov and N. N. Seregin, Metrologia 51, 54 (2014).
  * Flowers  _et al._ [1993] J. L. Flowers, B. W. Petley, and M. G. Richards, Metrologia 30, 75 (1993).
  * Aruev and Neronov [2012] N. N. Aruev and Y. Neronov, Tech. Phys. 57, 1579 (2012).
  * Jackowski  _et al._ [2010] K. Jackowski, M. Jaszuński, and M. Wilczek, The Journal of Physical Chemistry A 114, 2471 (2010).
  * Garbacz  _et al._ [2012] P. Garbacz, K. Jackowski, W. Makulski, and R. E. Wasylishen, J. Phys. Chem. A 116, 11896 (2012).
  * Tiesinga  _et al._ [2021] E. Tiesinga, P. J. Mohr, D. B. Newell, and B. N. Taylor, Rev. Mod. Phys. 93, 025010 (2021).
  * Schneider  _et al._ [2019] A. Schneider, A. Mooser, A. Rischka, K. Blaum, S. Ulmer, and J. Walz, Ann. Phys. 531, 1800485 (2019).
  * Shifrin  _et al._ [1997] V. Y. Shifrin, P. G. Park, C. G. Kim, V. N. Khorev, and C. H. Choi, IEEE Trans. Instrum. Meas. 46, 97 (1997).
  * Mooser  _et al._ [2018] A. Mooser, A. Rischka, A. Schneider, K. Blaum, S. Ulmer, and J. Walz, J. Phys.: Conf. Ser. 1138, 012004 (2018).
  * Mooser [2021] A. Mooser, private communication (2021).
  * Mooser  _et al._ [2014] A. Mooser, S. Ulmer, K. Blaum, K. Franke, H. Kracke, C. Leiterlitz, W. Quint, C. C. Rodegheri, C. Smorra, and J. Walz, Nature 509, 596 (2014).
  * Schneider  _et al._ [2017] G. Schneider, A. Mooser, M. Bohman, N. Schön, J. Harrington, T. Higuchi, H. Nagahama, S. Sellner, C. Smorra, K. Blaum, Y. Matsuda, W. Quint, J. Walz, and S. Ulmer, Science 358, 1081 (2017).
  * Smorra  _et al._ [2017] C. Smorra, S. Sellner, M. J. Borchert, J. A. Harrington, T. Higuchi, H. Nagahama, T. Tanaka, A. Mooser, G. Schneider, M. Bohman, K. Blaum, Y. M. C. Ospelkaus, W. Quint, J. Walz, Y. Yamazaki, and S. Ulmer, Nature 350, 371 (2017).
  * Lamb [1941] W. E. Lamb, Phys. Rev. 60, 817 (1941).
  * Ramsey [1950] N. F. Ramsey, Phys. Rev. 78, 699 (1950).
  * Caswell and Lepage [1986] W. E. Caswell and G. P. Lepage, Phys. Lett. B 167, 437 (1986).
  * Ivanov  _et al._ [2009] V. G. Ivanov, S. G. Karshenboim, and R. N. Lee, Phys. Rev. A 79, 012512 (2009).
  * Rudziński  _et al._ [2009] A. Rudziński, M. Puchalski, and K. Pachucki, J. Chem. Phys. 130, 244102 (2009).
  * Pachucki [2008] K. Pachucki, Phys. Rev. A 78, 012504 (2008).
  * Gimenez  _et al._ [2016] C. A. Gimenez, K. Kozioł, and G. A. Aucar, Phys. Rev. A 93, 032504 (2016).
  * Kutzelnigg [2012] W. Kutzelnigg, Chemical Physics 395, 16 (2012), recent Advances and Applications of Relativistic Quantum Chemistry.
  * Yerokhin  _et al._ [2011] V. A. Yerokhin, K. Pachucki, Z. Harman, and C. H. Keitel, Phys. Rev. Lett. 107, 043004 (2011).
  * Yerokhin  _et al._ [2012] V. A. Yerokhin, K. Pachucki, Z. Harman, and C. H. Keitel, Phys. Rev. A 85, 022512 (2012).
  * Pachucki [2004] K. Pachucki, Phys. Rev. A 69, 052502 (2004).
  * Itzykson and Zuber [2005] C. Itzykson and J.-B. Zuber,  _Quantum Field Theory_ (Dover Publications, New York, 2005).
  * Jentschura  _et al._ [2005] U. D. Jentschura, A. Czarnecki, and K. Pachucki, Phys. Rev. A 72, 062102 (2005).
  * Bethe and Salpeter [1977] H. A. Bethe and E. E. Salpeter,  _Quantum Mechanics of One- and Two-Electron Atoms_ (Plenum Publishing Corporation, New York, 1977).
  * Helgaker  _et al._ [1999] T. Helgaker, M. Jaszuński, and K. Ruud, Chemical Reviews 99, 293 (1999).
  * Pachucki [1998] K. Pachucki, J. Phys. B: At., Mol. Opt. Phys. 31, 5123 (1998).
  * Puchalski  _et al._ [2013] M. Puchalski, J. Komasa, and K. Pachucki, Phys. Rev. A 87, 030502(R) (2013).
  * Pachucki  _et al._ [2005] K. Pachucki, A. Czarnecki, U. D. Jentschura, and V. A. Yerokhin, Phys. Rev. A 72, 022108 (2005).
  * Aznabaev  _et al._ [2018] D. T. Aznabaev, A. K. Bekbaev, and V. I. Korobov, Phys. Rev. A 98, 012510 (2018).
  * Pachucki  _et al._ [2017] K. Pachucki, V. Patkóš, and V. A. Yerokhin, Phys. Rev. A 95, 062510 (2017).
  * Vaara and Pyykkö [2003] J. Vaara and P. Pyykkö, J. Chem. Phys. 118, 2973 (2003).
  * Kudo and Fukui [2005] K. Kudo and H. Fukui, J. Chem. Phys. 123, 114102 (2005).
  * Antušek  _et al._ [2007] A. Antušek, M. Jaszuński, and A. Rizzo, J. Chem. Phys. 126, 074303 (2007).
  * Seino and Hada [2010] J. Seino and M. Hada, J. Chem. Phys. 132, 174105 (2010).
  * Kupka  _et al._ [2013] T. Kupka, M. Stachów, L. Stobiński, and J. Kaminský, Magn. Reson. Chem. 51, 463 (2013).
  * [46] D.J. Wineland, W.M. Itano, R.S. van Dyck Jr., Adv. At. Mol. Phys. 19, 135 (1983).
  * [47] A.Antušek, P. Rodziewicz, D. Kȩdziera, A. Kaczmarek-Kȩdziera, M, Jaszuński, Chem. Phys. Lett. 588, 57 (2013).
  * [48] K. Pachucki and M. Puchalski, Optics Comm. 283, 641 (2010).