# Electron induced nanoscale nuclear spin relaxation probed by
hyperpolarization injection

William Beatrez  Department of Chemistry, University of California, Berkeley,
Berkeley, CA 94720, USA.  Arjun Pillai  Department of Chemistry, University of
California, Berkeley, Berkeley, CA 94720, USA.  Otto Janes  Department of
Chemistry, University of California, Berkeley, Berkeley, CA 94720, USA.
Dieter Suter  Fakultät Physik, Technische Universität Dortmund, D-44221
Dortmund, Germany.  Ashok Ajoy<EMAIL_ADDRESS>Department of Chemistry,
University of California, Berkeley, Berkeley, CA 94720, USA.  Chemical
Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720,
USA.

###### Abstract

We report on experiments that quantify the role of a central electronic spin
as a relaxation source for nuclear spins in its nanoscale environment. Our
strategy exploits hyperpolarization injection from the electron as a means to
controllably probe an increasing number of nuclear spins in the bath, and
subsequently interrogate them with high fidelity. Our experiments are focused
on a model system of a nitrogen vacancy (NV) center electronic spin surrounded
by several hundred ${}^{13}\textrm{C}$ nuclear spins. We observe that the
${}^{13}\textrm{C}$ transverse spin relaxation times vary significantly with
the extent of hyperpolarization injection, allowing the ability to measure the
influence of electron mediated relaxation extending over several nanometers.
These results suggest interesting new means to spatially discriminate nuclear
spins in a nanoscale environment, and have direct relevance to dynamic nuclear
polarization and quantum sensors and memories constructed from hyperpolarized
nuclei.

Introduction – Coupled electron-nuclear spin systems are highly relevant to
quantum information processing [1, 2], quantum sensing [3], and dynamic
nuclear polarization (DNP) [4]. Consider an electronic spin centrally located
(position $r{=}0$) in a bath of nuclear spins within a magnetic field $B_{0}$
(see Fig. 1A-B), coupled to them via hyperfine interactions. The electrons are
fast-relaxing (short $T_{1e}$), often on account of interactions with phonons
[5], and in turn they can serve as a relaxation source for the nuclear spins
[6, 7, 8] (schematically shown in Fig. 1C). Such electron induced nuclear
relaxation is an important consideration for applications in quantum
registers, memories, and sensors constructed out of nuclear spins [9, 10, 11],
and is therefore important to quantify. It also plays a key role in
determining rates of polarization transfer in DNP [12, 13, 14]. However,
probing such relaxation influences, particularly in a spatially defined
manner, is challenging. This is because experiments typically have very
restricted access to spins in the bath, i.e. there is a limited possibility of
spatially distinguishing the spins, other than a small shell where the nuclear
resonance frequencies can be significantly shifted from those of the bulk [8,
15, 16, 17, 18]. Previous nanoscale quantum sensing experiments have, for
instance, been limited to proximal central spin relaxation effects in small
($<$20) spin networks [19, 20].

In this paper, we report on experiments that study the effects of an
electronic spin on nuclei over wider length scales, spanning several
nanometers and involving several hundred nuclei. Our strategy (Fig. 1D)
exploits controllable hyperpolarized spin injection from the electron to the
nuclear bath. The hyperpolarization time $\tau$ is employed as a knob to tune
length scales in the bath being probed; the polarization is carried over
longer distances with increasing $\tau$. Simultaneously, the nuclei are
subject to a RF driving protocol that permits continuous tracking of their
magnetization over minute long periods with high signal-to-noise (SNR) [21].
The measurements reveal, surprisingly, that ensemble-averaged nuclear
lifetimes $T_{2}^{\prime}$ in the rotating frame increase dramatically with
increasing polarization time $\tau$. We demonstrate that this constitutes a
direct experimental signature of electron induced nuclear spin relaxation,
allowing us to quantify its influence extending over several nanometers.

Fig. 1: System and protocol.  (A-C) System consists of central NV electron
(red), ${}^{13}\textrm{C}$ nuclei (green), and P1 centers (blue) at distance
$r$. Dashed lines are inter-nuclear dipolar couplings. Profile denotes
$1/r^{6}$ NV-mediated relaxation effect. (A,C) Role of NV center as
polarization source and relaxation sink respectively, during different
experimental regimes. (D) Experiment schematic:(I)
NV${\rightarrow}{}^{13}\textrm{C}$ hyperpolarization for period $\tau$ at
36mT, (II) transport to high field and (III) ${}^{13}\textrm{C}$ readout for
time $t$ at 7T. (E) Hyperpolarization (I) involves MW chirps [22, 23]. (F)
Measurement (III) comprises a train of spin-locking $\vartheta$-pulses with
interrogation in $t_{\textrm{acq}}$ inter-pulse intervals [24]. Note:
Throughout the manuscript, time $\tau$ describes the hyperpolarization (spin-
lock) time and time $t$ describes the readout (detection) time.  Fig. 2:
Effect of increasing hyperpolarization time $\tau$.  (A) Spin-lock decays for
different hyperpolarization times $\tau$ corresponding to colorbar (regime I
in Fig. 1D). Single-shot data here (obtained with $t_{\textrm{acq}}{=}32\mu$s
and $\vartheta{\approx}\pi/2$ in Fig. 1F) is boxcar averaged over 97ms, and
normalized after truncation at $t{=}9.7$ms. Data for $t{=}36$s corresponds to
$\sim$369k pulses (see Fig. 3 for full data). Dashed line represents $1/e$
intercept. Measured $T_{2}^{\prime}$ increases ${\approx}6.42$-fold between
$\tau{=}2$s and $\tau{=}120s$. Inset (i): Zoom into representative 1s segment
(gray window in (A)). Dashed line is a fitted piecewise monoexponential. (B)
Instantaneous lifetime $T_{2}$ measured along decay curve, extracted from
slopes as in A(i). Color bar (see (A)) represents $\tau$. Dashed line shows
exemplary segment ending at 16s (gray window in (A)). (C) Instantaneous
$T_{2}$ plotted against $\tau$. Points show $T_{2}$ lifetimes from (B) for
exemplary 1s segments ending at labeled $t$ values (different markers). Lines
are spline fits. Fig. 3: Quantifying instantaneous decay rates.  (A) Full data
in Fig. 2 plotted on a logarithmic scale with respect to $\sqrt{t}$ for
different $\tau$ (color bar). Traces extending to $t{=}144$s have
${\sim}10^{6}$ points. Upper axis represents time $t$. Six equally sized (in
$\sqrt{t}$) segments (labeled (i)-(vi)) are selected along the line. (B) Zoom
into segment (i) (sampled every $\approx$5ms) for representative marked $\tau$
values (same colorbar as in (A)). Dashed lines are at equally spaced
$\sqrt{t}$ values to guide the eye. (C) Extracted time constants
$T_{2}^{\prime}$ from slopes of the corresponding segments in (A). Segment (i)
shows a steep increase in $T_{2}^{\prime}$ with $\tau$; subsequent segments
show progressively flatter profiles (gray arrow). Crossover of early and late
segments occurs at $\tau{\approx}32$s (dashed line). (D) Variations in short
and long time decay behavior with $\tau$. At long $t$ (dark shaded region),
decay closely follows $\alpha{=}1/2$ (dashed lines, extracted from segment
(vi)), while short times $t$ show a transition from convex to concave
behavior. Transition occurs at $\tau{\approx}32s$ (similar to Fig. 3C). Upper
axis represents time $t$. Line colors correspond to colorbar in panel (A).
Fig. 4:  Simulations of polarization evolution in time and space. (A) Spatial
polarization build-up with $\tau$ (colorbar) plotting polarization
$\mathcal{P}(r,\tau)$ in a shell at distance $r$ from the central NV center
(during regime I in Fig. 1D). (B) Polarization distribution during readout
period $t$ (colorbar) in regime III of Fig. 1D. Here, we set the initial
condition to be the $\tau{=}120$s polarization distribution in panel (A). (C)
Simulated NMR signals following panels (A-B) for different hyperpolarization
times $\tau$ similar to Fig. 2A. Data reveals $T_{2}^{\prime}$ increase with
$\tau$, in agreement with the experimental data in Fig. 2 (see (A) for
colorbar). (D) Stretched exponential decay dynamics. Analogous to Fig. 3A,
data from (C) is plotted on a logarithmic scale with respect to $\sqrt{t}$.
Stretched exponential behavior qualitatively matches data in Fig. 3A. Segments
labeled (i)-(vi) are used to extract $T_{2}^{\prime}$ decay constants. (E)
Extracted $T_{2}^{\prime}$ decay constants for representative segments
(i)-(vi) show qualitative agreement with Fig. 3C. Here, dark lines are
simulation and light lines are experimental data. Dashed line marks crossover
at $\tau{\sim}32$s.

System and Protocol – Experiments here are carried out in diamond as a model
system, with central NV [25, 26] electronic spins surrounded by
${}^{13}\textrm{C}$ nuclei (Fig. 1A-B). At ${\approx}$1ppm, the NV centers are
separated by ${\approx}12$nm, and ${}^{13}\textrm{C}$ nuclei appear with
lattice density ${\approx}$0.92 spins/nm3, yielding a relative
NV:${}^{13}\textrm{C}$ ratio ${\sim}1{:}10^{4}$ [27, 28]. In addition, the
lattice hosts randomly positioned paramagnetic impurities (P1 centers) at a
concentration $>$20ppm, which can also serve as relaxation sources (Fig.
1A-B).

Experiments are conducted in three regimes (indicated in Fig. 1D): (I)
optically induced NV${\rightarrow}{}^{13}\textrm{C}$ spin injection for period
$\tau$ at low-field ($B_{\textrm{pol}}{=}36$mT), (II) rapid adiabatic
transport to high field, and (III) ${}^{13}\textrm{C}$ interrogation at
$B_{0}{=}7$T. For hyperpolarization (Fig. 1E), we exploit a mechanism
previously described in Ref. [22, 23, 29]. Distant spins are polarized by spin
diffusion driven by the internuclear dipolar Hamiltonian [18, 30, 31],
$\mathcal{H}_{\textrm{dd}}{=}\sum_{j<k}d_{jk}^{\textrm{CC}}(3I_{jz}I_{kz}-\vec{I_{j}}\cdot\vec{I_{k}})$,
where $I_{j}$ refer to $j^{\textrm{th}}$ spin-1/2 operators, and coupling
strengths $d_{jk}^{\textrm{CC}}{\propto}\gamma_{n}^{2}/r^{3}$, with
magnetogyric ratio $\gamma_{n}{=}$10.7MHz/T, and $r$ being the inter-nuclear
distance. Bulk averaged ${}^{13}\textrm{C}$ hyperpolarization under maximal
conditions is $\varepsilon{\approx}0.3\%$.

Sample transport (regime II) occurs in $t_{s}{\approx}1\textrm{s}{\ll}T_{1n}$
in a manner adiabatic with respect to the instantaneous ${}^{13}\textrm{C}$
Larmor frequency, thereby preserving the hyperpolarization generated [32].
Conversely, since $t_{s}{\gg}T_{1e}{\sim}1$ms [33, 5], the NV center rapidly
loses hyperpolarization and ultimately returns to thermal polarization levels
${\approx}3$% in regime III (Fig. 1F). Subsequently, the NV center
predominantly serves as a point relaxation source for the nuclear bath (Fig.
1C).

Ensemble ${}^{13}\textrm{C}$ readout (regime III) employs a protocol described
in Refs. [34, 24]. ${}^{13}\textrm{C}$ spins are prepared along the transverse
axis $\hat{\textbf{x}}$ ($\rho_{I}{=}\varepsilon\sum_{j}I_{jx}$) on the Bloch
sphere, and a train of spin-locking $\vartheta$-pulses are applied [35, 36].
${}^{13}\textrm{C}$ nuclei are interrogated in windows between the pulses,
allowing their dynamics to be continuously tracked with high SNR [37]. Signal
obtained corresponds to measuring the survival probability in the
$\hat{\textbf{x}}{\text{-}}\hat{\textbf{y}}$ plane. The sequence operation can
be described by the unitary $U(t)$. For sufficiently rapid pulsing duty cycle,
$U(t){\approx}\exp(i{\mathcal{H}}^{(0)}_{F}t)$, such that the inter-nuclear
Hamiltonian is engineered to leading order in the Magnus expansion to
${\mathcal{H}}^{(0)}_{F}{\approx}\sum_{j<k}d_{jk}^{\textrm{CC}}\left(\frac{3}{2}\mathcal{H}_{\textrm{ff}}-\vec{I_{j}}\cdot\vec{I_{k}}\right)$,
where $\mathcal{H}_{\textrm{ff}}{=}I_{jz}I_{kz}+I_{jy}I_{ky}$ is a flip-flop
Hamiltonian [24]. Since $[\rho_{I}{,}{\mathcal{H}}^{(0)}_{F}]{=}0$ commutes
with the initial state, dipolar evolution is suppressed to leading order. As a
result, ${}^{13}\textrm{C}$ free induction decay lifetimes
$T_{2}^{*}{\sim}1.5$ms (in the absence of spin-lock pulses), are significantly
extended, here to $T_{2}^{\prime}{\gtrsim}$65.5s (see Fig. 2A).

Results — Fig. 2A describes our primary experimental result (see movie at Ref.
[38]), showing ${}^{13}\textrm{C}$ NMR signal measured employing differing
hyperpolarization periods $\tau$ (varied every one second from $\tau{=}$1s to
$\tau{=}120$s, see colorbar). The signals here are normalized to their values
at 9.7ms. Each 36s trace consists of ${\sim}369,000$ pulses (full data in Fig.
3 consists of ${\sim}{10^{6}}$ pulses), and the ${}^{13}\textrm{C}$ nuclei are
interrogated after every pulse (Fig. 1F). Surprisingly, we observe that the
signals decay more slowly with increasing hyperpolarization period $\tau$.
Normalization of the signal profiles allows the ability to unravel the
relative changes in the decay time constant $T_{2}^{\prime}$, estimated from
the $1/e$ intercept (dashed line in Fig. 2A). For example, comparing
$\tau{=}2$s and $\tau{=}$120s in Fig. 2A, we observe a $T_{2}^{\prime}$
increase from 10.2s to 65.5s.

The $1/e$-intercept (dashed line in Fig. 2A) is blind to the instantaneous
change of the decay profile and therefore provides only limited information.
To more clearly observe the decay dynamics, we divide the curves in Fig. 2A
into 36 segments of width $\Delta{t}{=}1$s, with one such segment shown in the
gray window in Fig. 2A. We fit the decay in each segment to a monoexponential
(as in inset (i)) and extract the instantaneous time constants $T_{2}(t)$,
which are plotted in Fig. 2B. The notation $T_{2}$ (as opposed to
$T_{2}^{\prime}$) is used to emphasize that these are monoexponential
constants. Colors here represent $\tau$ with the same color bar as in Fig. 2A.
For each trace, the signals decay markedly slower with increasing time $t$
(see also Fig. 3). Ultimately, the $T_{2}$ times are remarkably long
($T_{2}{\approx}$100s) at large ${\tau}$. Increasing polarization time $\tau$
makes the overall decay slower for any selected segment. Indeed, for the
segment ending at $t{=}16$s (vertical dashed line in Fig. 2B), the $T_{2}$
value is increased by $1.66$-fold. To now emphasize the relative change in the
$T_{2}$ values for different segments, Fig. 2C shows the $T_{2}$ lifetimes
plotted against $\tau$, where data corresponding to each segment in Fig. 2A-B
forms the lines (values denote segments ending at $t$). Fig. 2C makes clear
that the instantaneous $T_{2}$ increases for each segment with increasing $t$,
and within each segment with hyperpolarization time $\tau$. The maximum
relative change occurs for short $\tau$ and at early $t$.

A clearer view of data in Fig. 2 can be obtained by noting that the decays
approximately follow a stretched exponential of the form
$\exp\left[-(t/T_{2}^{\prime})^{\alpha}\right]$ with $\alpha{\approx}1/2$ (see
Fig. 3 and movies at Refs. [39, 40]). Fig. 3A shows the ${}^{13}\textrm{C}$
signals in Fig. 2 plotted on a logarithmic scale with respect to $\sqrt{t}$;
the colorbar represents increasing hyperpolarization time $\tau$. The signals
(gray lines) then manifest as approximately straight lines, demonstrating
stretched exponential character with $\alpha{\approx}1/2$. $T_{2}^{\prime}$
lifetimes can now be extracted conveniently from the instantaneous slopes $s$,
as $T_{2}^{\prime}{=}1/s^{2}$. High SNR, along with rapid data collection
rates (allowing ${\sim}10^{6}$ points per trace), allow the unique ability to
continuously observe the stretched exponential dynamics. We now focus
attention to six segments along the decay curves (labeled (i)-(vi) in Fig.
3A). Complementary to Fig. 2A, data reveals that the $T_{2}^{\prime}$ values
at short time $t$ vary considerably with increasing $\tau$, as evident in
segment (i). In contrast, the $T_{2}^{\prime}$ values at long time $t$ are
observed to be independent of $\tau$ (evident from the approximately parallel
relaxation profiles in segment (vi)). Segment (i) is zoomed in Fig. 3B on the
same logarithmic scale for representative values of $\tau$. The panel makes
clear the large relative change in $T_{2}^{\prime}$ with $\tau$, most evident
in the difference between traces corresponding to $\tau{=}$1s and
$\tau{=}$120s.

To quantify this change with greater clarity, Fig. 3C elucidates the extracted
$T_{2}^{\prime}$ values for the six segments, with a special emphasis on
segments (i) (0-1s) and (vi) (121-144s) (highlighted traces). Data reveals
that the $T_{2}^{\prime(i)}$ grows significantly with increasing
hyperpolarization time $\tau$, while $T_{2}^{\prime(vi)}$ is flat and almost
independent of $\tau$. Indeed, the $\tau$ dependence of $T_{2}^{\prime}$
becomes “flatter” with increasing segment number. Interestingly also, the
extracted $T_{2}$ values for the different segments display a crossing point
at $\tau{\approx}$32s (dashed line in Fig. 3C). Finally, Fig. 3D offers
insight into variations from exact $\alpha{=}1/2$ stretched exponential
behavior with changing $\tau$. We observe that at long times $t$ (dark shaded
region), the behavior follows a universal $\alpha{=}1/2$ dependence (dashed
lines). However, in the short time region ($t{\lesssim}$36s), there is a
deviation from this, exhibited as a transition from convex to concave behavior
around $\alpha{=}1/2$. Cross-over occurs at $\tau{\approx}$32s, similar to
Fig. 3C.

While similar stretched exponential decay dynamics have been observed before
in some restrictive model systems [41, 42, 43], the relatively low number of
measurement points in these reports allowed only a limited ability to
characterize them. By contrast, the ${\sim}10^{6}$ points per trace in Fig.
3A, along with the high SNR, provide a remarkably clear view into the dynamics
at short and long $t$.

Theory – To describe the experimental observations, we construct a simple,
semi-quantitative model for the nuclear polarization $p(r,t)$ at coordinate
$r$ and time $t$. We assume centrosymmetry, a good approximation given the
ensemble average in our experiments. We then model the dynamics using the
differential equation,

$\frac{\partial}{\partial
t}p(r,t)=\frac{P_{0}}{r^{6}}-\frac{\kappa_{0}}{r^{6}}p(r,t)-\frac{1}{T_{1}}p(r,t)+D\nabla^{2}p(r,t)$
(1)

where $P_{0}$ denotes the rate of hyperpolarization injection and $\kappa_{0}$
is the strength of spatially dependent relaxation due to the central NV
center. The $1/r^{6}$ scaling of the $P_{0}$ term in Eq. (1) does not yield
qualitatively different behavior compared to a $1/r^{3}$ scaling but makes the
equation better conditioned near $r{=}0$. In any case, the $r{\rightarrow}0$
spatial regime and associated nuclear spin frozen core are outside the ambit
of our experimental observations because of the relatively small
${}^{13}\textrm{C}$ detection bandwidth (${\approx}$32kHz) in our
measurements. In contrast to the NV center (at $r{=}0$), we assume that the
relatively dense P1 centers (${>}$20ppm) serve as contributors to background
relaxation of the ${}^{13}\textrm{C}$ nuclei independent of their position;
this is captured by the $T_{1}$ term in Eq. (1). Finally, the last term
denotes spin diffusion which we assume to be Fickian and Gaussian with a
single constant $D$ at all locations. This is a good approximation given the
large number of ${}^{13}\textrm{C}$ nuclei participating in the measured
signal around every NV center [18]. For instance, a sphere of injected
polarization with radius 4nm contains ${\approx}247$ nuclear spins.

We solve Eq. (1) separately in the three regimes of the experiment (Fig. 1D),
with the solution for one regime setting the initial conditions for the next.
Obtaining a precise estimate of the parameters in Eq. (1) is challenging and
outside the scope of this paper; our goals instead are to obtain qualitative
agreement with the experimental observations in Fig. 2-Fig. 3. We therefore
make some simplifying assumptions to the parameters in Eq. (1). In regime I,
we assume $P_{0}/\kappa_{0}{=}1/100$ for simplicity ($P_{0}{=}0$ in regime
III). In regime I and III, we employ $T_{1}{=}q^{-1}T_{1,\textrm{LF}}$ and
$T_{1}{=}q^{-1}T_{1,\textrm{HF}}$ respectively, where
$T_{1,\textrm{LF}}{=}283$s and $T_{1,\textrm{HF}}{=}1520$s are the measured
low and high-field bulk ${}^{13}\textrm{C}$ lifetimes respectively at
$\tau{=}60$s [37], and $q$ is a scaling factor that we employ as a free
parameter in the fits. For Fig. 4, we find good agreement with $q{=}6.75$. The
latter assumption can be rationalized by the fact that (1) the individual
contributions from P1 centers to $T_{1}$ relaxation are hard to separate in
the bulk $T_{1}$ measurements; (2) longitudinal and transverse relaxation
rates in regime III are expected to be proportional, but are in general,
unequal, making the measured $T_{1,\textrm{HF}}$ only approximately reflective
of $T_{1}$ in Eq. (1); and (3) $T_{1,\textrm{LF}}$ is measured under dark
conditions at a fixed $\tau{=}60$s, but the corresponding $T_{1}$ in regime I
is measured under optical illumination (making measurement more challenging).

With these assumptions, Fig. 4A-B shows the simulated system dynamics, where
we plot the polarization contained in a shell at radius $r$,
$\mathcal{P}{=}4\pi r^{2}p$. Fig. 4A first shows $\mathcal{P}(r,\tau)$ with
increasing hyperpolarization time $\tau$ (see colorbar) in regime I, assuming
a 1s shuttling period in regime II. Simulation parameters are set to obtain
good qualitative agreement with experiment (here $D{=}0.0135$nm2/s and
$P_{0}{=}6.75$s-1). With increasing $\tau$, spin diffusion leads to a spread
of polarization; this is evidenced by the movement of the ”centroid” of
$\mathcal{P}(r,\tau)$ with increasing $\tau$ towards greater $r$ in Fig. 4A.
Notably however, the continuous replenishment of polarization from the NV in
regime I makes the distribution of $\mathcal{P}(r,\tau)$ skewed towards the
left. In a complementary manner, Fig. 4B elucidates the distribution
$\mathcal{P}(r,t)$ during the readout period in regime III, assuming one
starts with the distribution obtained with $\tau{=}120$s (right-most trace in
Fig. 4A). The strong relaxing effect of the NV center yields the polarization
“hole” close to $r{=}0$, and manifests as the steep wall of growing
polarization in Fig. 4B. Additionally, the centroid of $\mathcal{P}(r,t)$
moves towards larger $r$ and homogenizes as $t$ increases; the shift with $r$
is apparently larger here because there is no polarization replenishment and
relaxation is slower for larger $r$.

Fig. 4C displays the net polarization $\int\mathcal{P}(r,t)dr$ following the
trajectory of Fig. 4B as interrogated during the readout period, but
considering here the result for different hyperpolarization times $\tau$ (see
colorbar in Fig. 4A). Normalizing the traces then yields an apparent
polarization relaxation that becomes slower with increasing $\tau$, closely
matching the experimental results in Fig. 2A. From Fig. 4A-B, we identify that
the observations in Fig. 2A arise because the shifting centroid of
$\mathcal{P}$ (see Fig. 4A) makes electron mediated relaxation less effective
with increasing $\tau$. Interestingly, we find that the decays in Fig. 4A also
follow an identical $\alpha{=}1/2$ stretched exponential as in the experiment.
This is evidenced by replotting the data in a logarithmic scale with respect
to $\sqrt{t}$ in Fig. 4D, where we observe a behavior similar to the
experiments in Fig. 3A: the stretched exponential dynamics manifest as
straight lines here. We hypothesize that the stretched exponentials here
result from the inhomogeneity of the polarization distribution in Fig. 4A,
which is subject to a collection of relaxation rates based on spatial
proximity to the NV. At long readout periods $t$, the parallel lines yield a
universal relaxation profile that is independent of $\tau$; Fig. 4B allows us
to recognize that this is because the polarization spreading far from the NV
center is predominantly influenced by the background $1/T_{1}$ relaxation.
Indeed, as shown in Fig. 4E, taking segments (i)-(vi) along the decay curves
in Fig. 4D, and extracting their relaxation times $T_{2}^{\prime}$, we observe
a progressive flattening of the $T_{2}^{\prime}$ values with increasing
$\tau$, in good agreement with experimental data in Fig. 3C (overlapped here).
Cross-over occurs at $\tau{\sim}32$s similar to Fig. 3C.

Outlook – These results open many interesting future avenues in a wide class
of central spin systems. First, as Fig. 3A and Fig. 4 demonstrate, spin-lock
control allows for peering into nuclear polarization localized at different
positions $r$ with respect to the central electron. Moving through time $t$,
the decay curve, as seen in Fig. 2A, effectively corresponds to shifting the
centroid location of the sensitive region being probed in the lattice (Fig.
4B). This suggests a (nonlinear) means to map from $t$ to an effective $r$
coordinate, and can allow interesting new ways to discriminate nuclear spins
in the electronic environment. It can enable visualization of the transport of
polarization from the central electron [44, 45, 46, 47], and the ability to
rewritably engineer spin texture in the nuclear bath [48, 49]. For instance,
successive injection of positive and negative polarization from the NV center
can result in non-equilibrium “shells” of nuclear polarization whose dynamics
can be frozen for minute-long periods by decoupling. Ultimately, this
engenders opportunities for exploiting mesoscale nuclear spin baths for
quantum information science (QIS) [50] and sensing [51, 52, 53, 54].

Our results illustrate the ability of nuclear spins to probe electronic
relaxation processes, therefore providing a means to view the phonon density
of states that dominate these relaxation mechanisms [5, 55]. This is relevant
to the emerging field of QIS with molecular systems [56], where there is an
important need to develop design rules for electronic lifetimes based on
molecular vibronic properties [57]. It also has important implications for
DNP, wherein the concentration and identity of electron spins can vastly
affect nuclear polarization levels [58]. Concurrently, this suggests exciting
potential for leveraging long nuclear $T_{2}^{\prime}$ lifetimes in systems
capable of “turning off” the electron spin after hyperpolarization injection
[59], for instance in molecules hosting triplet photoexcitable electrons [60,
61, 62].

Conclusion – In summary, we have experimentally measured the influence of a
central electronic spin in relaxing the surrounding nuclear spin bath. We
exploited controlled hyperpolarization injection as a knob to control length
scales over which bath spins are probed, extending over several nanometers and
hundreds of spins. Our work informs on a new means to probe nanoscale spin
environments, and portends applications in quantum memories and sensors
constructed out of hyperpolarized nuclei.

We gratefully acknowledge discussions with C. Ramanathan, C. Meriles and J.
Reimer. This work was funded by ONR (N00014-20-1-2806) and DOE STTR (DE-
SC0022441).

## References

  * Prokof’ev and Stamp [2000] N. V. Prokof’ev and P. C. E. Stamp, Theory of the spin bath, Reports on Progress in Physics 63, 669 (2000).
  * Ladd  _et al._ [2010] T. Ladd, F. Jelezko, R. Laflamme, Y. Nakamura, C. Monroe, _et al._ , Quantum computers, Nature 464, 45 (2010).
  * Degen  _et al._ [2017] C. L. Degen, F. Reinhard, and P. Cappellaro, Quantum sensing, Reviews of modern physics 89, 035002 (2017).
  * Abragam and Goldman [1978] A. Abragam and M. Goldman, Principles of dynamic nuclear polarization, Reports on Progress in Physics 41, 395 (1978).
  * Jarmola  _et al._ [2012] A. Jarmola, V. Acosta, K. Jensen, S. Chemerisov, and D. Budker, Temperature-and magnetic-field-dependent longitudinal spin relaxation in nitrogen-vacancy ensembles in diamond, Physical review letters 108, 197601 (2012).
  * Witzel and Sarma [2006] W. M. Witzel and S. D. Sarma, Quantum theory for electron spin decoherence induced by nuclear spin dynamics in semiconductor quantum computer architectures: Spectral diffusion of localized electron spins in the nuclear solid-state environment, Phys. Rev. B 74, 035322 (2006).
  * Bloembergen  _et al._ [1947] N. Bloembergen, E. Purcell, and R. Pound, Nuclear magnetic relaxation, Nature 160, 475 (1947).
  * Blumberg [1960] W. Blumberg, Nuclear spin-lattice relaxation caused by paramagnetic impurities, Physical Review 119, 79 (1960).
  * Hanson  _et al._ [2008] R. Hanson, V. V. Dobrovitski, A. E. Feiguin, O. Gywat, and D. D. Awschalom, Coherent dynamics of a single spin interacting with an adjustable spin bath, Science 320, 352 (2008).
  * Reiserer  _et al._ [2016] A. Reiserer, N. Kalb, M. S. Blok, K. J. van Bemmelen, T. H. Taminiau, R. Hanson, D. J. Twitchen, and M. Markham, Robust quantum-network memory using decoherence-protected subspaces of nuclear spins, Physical Review X 6, 021040 (2016).
  * Morton  _et al._ [2008] J. J. L. Morton, A. M. Tyryshkin, R. M. Brown, S. Shankar, B. W. Lovett, A. Ardavan, T. Schenkel, E. E. Haller, J. W. Ager, and S. A. Lyon, Solid-state quantum memory using the 31p nuclear spin, Nature 455, 1085 (2008).
  * Siaw  _et al._ [2014] T. A. Siaw, M. Fehr, A. Lund, A. Latimer, S. A. Walker, D. T. Edwards, and S.-I. Han, Effect of electron spin dynamics on solid-state dynamic nuclear polarization performance, Physical Chemistry Chemical Physics 16, 18694 (2014).
  * Saliba  _et al._ [2017] E. P. Saliba, E. L. Sesti, F. J. Scott, B. J. Albert, E. J. Choi, N. Alaniva, C. Gao, and A. B. Barnes, Electron decoupling with dynamic nuclear polarization in rotating solids, Journal of the American Chemical Society 139, 6310 (2017).
  * Sarkar  _et al._ [2021] A. Sarkar, B. Blankenship, E. Druga, A. Pillai, R. Nirodi, S. Singh, A. Oddo, P. Reshetikhin, and A. Ajoy, Rapidly enhanced spin polarization injection in an optically pumped spin ratchet, arXiv preprint arXiv:2112.07223 (2021).
  * Tse and Hartmann [1968] D. Tse and S. Hartmann, Nuclear spin-lattice relaxation via paramagnetic centers without spin diffusion, Physical Review Letters 21, 511 (1968).
  * Horvitz [1971] E. P. Horvitz, Nuclear spin diffusion induced by paramagnetic impurities in nonconducting solids, Physical Review B 3, 2868 (1971).
  * Stern  _et al._ [2021] Q. Stern, S. F. Cousin, F. Mentink-Vigier, A. C. Pinon, S. J. Elliott, O. Cala, and S. Jannin, Direct observation of hyperpolarization breaking through the spin diffusion barrier, Science Advances 7, eabf5735 (2021).
  * Ramanathan [2008] C. Ramanathan, Dynamic nuclear polarization and spin diffusion in nonconducting solids, Applied Magnetic Resonance 34, 409 (2008).
  * Dréau  _et al._ [2014] A. Dréau, P. Jamonneau, O. Gazzano, S. Kosen, J.-F. Roch, J. Maze, and V. Jacques, Probing the dynamics of a nuclear spin bath in diamond through time-resolved central spin magnetometry, Physical review letters 113, 137601 (2014).
  * Abobeih  _et al._ [2019] M. Abobeih, J. Randall, C. Bradley, H. Bartling, M. Bakker, M. Degen, M. Markham, D. Twitchen, and T. Taminiau, Atomic-scale imaging of a 27-nuclear-spin cluster using a quantum sensor, Nature 576, 411 (2019).
  * Sahin  _et al._ [2022] O. Sahin, H. A. Asadi, P. Schindler, A. Pillai, E. Sanchez, M. Elo, M. McAllister, E. Druga, C. Fleckenstein, M. Bukov, _et al._ , Continuously tracked, stable, large excursion trajectories of dipolar coupled nuclear spins, arXiv preprint arXiv:2206.14945 (2022).
  * Ajoy  _et al._ [2018a] A. Ajoy, K. Liu, R. Nazaryan, X. Lv, P. R. Zangara, B. Safvati, G. Wang, D. Arnold, G. Li, A. Lin, _et al._ , Orientation-independent room temperature optical 13c hyperpolarization in powdered diamond, Sci. Adv. 4, eaar5492 (2018a).
  * Ajoy  _et al._ [2018b] A. Ajoy, R. Nazaryan, K. Liu, X. Lv, B. Safvati, G. Wang, E. Druga, J. Reimer, D. Suter, C. Ramanathan, _et al._ , Enhanced dynamic nuclear polarization via swept microwave frequency combs, Proceedings of the National Academy of Sciences 115, 10576 (2018b).
  * Beatrez  _et al._ [2021] W. Beatrez, O. Janes, A. Akkiraju, A. Pillai, A. Oddo, P. Reshetikhin, E. Druga, M. McAllister, M. Elo, B. Gilbert, D. Suter, and A. Ajoy, Floquet prethermalization with lifetime exceeding 90 s in a bulk hyperpolarized solid, Phys. Rev. Lett. 127, 170603 (2021).
  * Jelezko and Wrachtrup [2006] F. Jelezko and J. Wrachtrup, Single defect centres in diamond: A review, Physica Status Solidi (A) 203, 3207 (2006).
  * Manson  _et al._ [2006] N. B. Manson, J. P. Harrison, and M. J. Sellars, Nitrogen-vacancy center in diamond: Model of the electronic structure and associated dynamics, Phys. Rev. B 74, 104303 (2006).
  * Ajoy  _et al._ [2019a] A. Ajoy, B. Safvati, R. Nazaryan, J. Oon, B. Han, P. Raghavan, R. Nirodi, A. Aguilar, K. Liu, X. Cai, _et al._ , Hyperpolarized relaxometry based nuclear t 1 noise spectroscopy in diamond, Nature communications 10, 1 (2019a).
  * Van Wyk  _et al._ [1997] J. Van Wyk, E. Reynhardt, G. High, and I. Kiflawi, The dependences of esr line widths and spin-spin relaxation times of single nitrogen defects on the concentration of nitrogen defects in diamond, Journal of Physics D: Applied Physics 30, 1790 (1997).
  * Ajoy  _et al._ [2020a] A. Ajoy, R. Nazaryan, E. Druga, K. Liu, A. Aguilar, B. Han, M. Gierth, J. T. Oon, B. Safvati, R. Tsang, _et al._ , Room temperature “optical nanodiamond hyperpolarizer”: Physics, design, and operation, Review of Scientific Instruments 91, 023106 (2020a).
  * Karabanov  _et al._ [2015] A. Karabanov, D. Wiśniewski, I. Lesanovsky, and W. Köckenberger, Dynamic nuclear polarization as kinetically constrained diffusion, Physical review letters 115, 020404 (2015).
  * Hovav  _et al._ [2010] Y. Hovav, A. Feintuch, and S. Vega, Theoretical aspects of dynamic nuclear polarization in the solid state–the solid effect, Journal of Magnetic Resonance 207, 176 (2010).
  * Ajoy  _et al._ [2019b] A. Ajoy, X. Lv, E. Druga, K. Liu, B. Safvati, A. Morabe, M. Fenton, R. Nazaryan, S. Patel, T. F. Sjolander, J. A. Reimer, D. Sakellariou, C. A. Meriles, and A. Pines, Wide dynamic range magnetic field cycler: Harnessing quantum control at low and high fields, Review of Scientific Instruments 90, 013112 (2019b), https://doi.org/10.1063/1.5064685 .
  * Popa  _et al._ [2004] I. Popa, T. Gaebel, M. Domhan, C. Wittmann, F. Jelezko, and J. Wrachtrup, Energy levels and decoherence properties of single electron and nuclear spins in a defect center in diamond, Phys. Rev. B 70, 201203 (2004).
  * Ajoy  _et al._ [2020b] A. Ajoy, R. Nirodi, P. Reshetikhin, E. Druga, A. Akkiraju, M. McAllister, G. Maineri, S. Le, A. Lin, A. M. Souza, C. A. Meriles, B. Gilbert, D. Suter, J. A. Reimer, and A. Pines, Dynamical decoupling in interacting systems: applications to signal-enhanced hyperpolarized readout, arXiv preprint arXiv:2008.08323 (2020b).
  * Rhim  _et al._ [1976] W.-K. Rhim, D. Burum, and D. Elleman, Multiple-pulse spin locking in dipolar solids, Physical Review Letters 37, 1764 (1976).
  * Rhim  _et al._ [1973] W.-K. Rhim, D. D. Elleman, and R. W. Vaughan, Analysis of multiple pulse nmr in solids, The Journal of Chemical Physics 59, 3740 (1973).
  * Beatrez  _et al._ [2022] W. Beatrez, A. Pillai, O. Janes, D. Suter, and A. Ajoy (2022), see Supplemental Material.
  * SD_[2021a] Video showing full dataset of Fig. 2A: https://www.youtube.com/watch?v=GEdVRV2cXH4 (2021a).
  * SD_[2021b] Video showing full dataset of Fig. 3A: https://www.youtube.com/watch?v=b8Be0R2etZQ (2021b).
  * SD_[2021c] Video showing unnormalized dataset of Fig. 3A: https://www.youtube.com/watch?v=5tirO0yE4zg (2021c).
  * Lock  _et al._ [1992] H. Lock, G. E. Maciel, and C. E. Johnson, Natural-abundance 13c dynamic nuclear polarization experiments on chemical vapor deposited diamond film, Journal of materials research 7, 2791 (1992).
  * Hartman  _et al._ [1994] J. S. Hartman, A. Narayanan, and Y. Wang, Spin-lattice relaxation in the 6h polytype of silicon carbide, Journal of the American Chemical Society 116, 4019 (1994).
  * Henrichs and Linder [1984] P. M. Henrichs and M. Linder, Carbon-13 spin diffusion in the determination of intermolecular structure in solids, Journal of Magnetic Resonance (1969) 58, 458 (1984).
  * Zhang and Cory [1998] W. Zhang and D. Cory, First direct measurement of the spin diffusion rate in a homogenous solid, Physical review letters 80, 1324 (1998).
  * Burgarth and Ajoy [2017] D. Burgarth and A. Ajoy, Evolution-free hamiltonian parameter estimation through zeeman markers, Physical Review Letters 119, 030402 (2017).
  * Ajoy  _et al._ [2012] A. Ajoy, R. K. Rao, A. Kumar, and P. Rungta, Algorithmic approach to simulate hamiltonian dynamics and an nmr simulation of quantum state transfer, Physical Review A 85, 030303 (2012).
  * Boutis  _et al._ [2004] G. Boutis, D. Greenbaum, H. Cho, D. Cory, and C. Ramanathan, Spin diffusion of correlated two-spin states in a dielectric crystal, Physical review letters 92, 137201 (2004).
  * King  _et al._ [2012] J. P. King, Y. Li, C. A. Meriles, and J. A. Reimer, Optically rewritable patterns of nuclear magnetization in gallium arsenide, Nature communications 3, 1 (2012).
  * Michal and Tycko [1999] C. A. Michal and R. Tycko, Stray-field nmr imaging and wavelength dependence of optically pumped nuclear spin polarization in inp, Physical Review B 60, 8672 (1999).
  * Foletti  _et al._ [2009] S. Foletti, H. Bluhm, D. Mahalu, V. Umansky, and A. Yacoby, Universal quantum control of two-electron spin quantum bits using dynamic nuclear polarization, Nature Physics 5, 903 (2009).
  * Sahin  _et al._ [2021] O. Sahin, E. d. L. Sanchez, S. Conti, A. Akkiraju, P. Reshetikhin, E. Druga, A. Aggarwal, B. Gilbert, S. Bhave, and A. Ajoy, High-field magnetometry with hyperpolarized nuclear spins, arXiv preprint arXiv:2112.11612 (2021).
  * Ajoy and Cappellaro [2012] A. Ajoy and P. Cappellaro, Stable three-axis nuclear-spin gyroscope in diamond, Physical Review A 86, 062104 (2012).
  * Jaskula  _et al._ [2019] J.-C. Jaskula, K. Saha, A. Ajoy, D. J. Twitchen, M. Markham, and P. Cappellaro, Cross-sensor feedback stabilization of an emulated quantum spin gyroscope, Physical Review Applied 11, 054010 (2019).
  * Ledbetter  _et al._ [2012] M. Ledbetter, K. Jensen, R. Fischer, A. Jarmola, and D. Budker, Gyroscopes based on nitrogen-vacancy centers in diamond, Physical Review A 86, 052116 (2012).
  * Takahashi  _et al._ [2008] S. Takahashi, R. Hanson, J. Van Tol, M. S. Sherwin, and D. D. Awschalom, Quenching spin decoherence in diamond through spin bath polarization, Physical review letters 101, 047601 (2008).
  * Bayliss  _et al._ [2020] S. Bayliss, D. Laorenza, P. Mintun, B. Kovos, D. Freedman, and D. Awschalom, Optically addressable molecular spins for quantum information processing, Science 370, 1309 (2020).
  * He  _et al._ [2020] C. He, P. R. Christensen, T. J. Seguin, E. A. Dailing, B. M. Wood, R. K. Walde, K. A. Persson, T. P. Russell, and B. A. Helms, Conformational entropy as a means to control the behavior of poly (diketoenamine) vitrimers in and out of equilibrium, Angewandte Chemie International Edition 59, 735 (2020).
  * Lange  _et al._ [2012] S. Lange, A. H. Linden, Ü. Akbey, W. T. Franks, N. M. Loening, B.-J. van Rossum, and H. Oschkinat, The effect of biradical concentration on the performance of dnp-mas-nmr, Journal of magnetic resonance 216, 209 (2012).
  * Capozzi  _et al._ [2017] A. Capozzi, T. Cheng, G. Boero, C. Roussel, and A. Comment, Thermal annihilation of photo-induced radicals following dynamic nuclear polarization to produce transportable frozen hyperpolarized 13c-substrates, Nature communications 8, 1 (2017).
  * Henstra and Wenckebach [2014] A. Henstra and W. T. Wenckebach, Dynamic nuclear polarisation via the integrated solid effect i: theory, Molecular Physics 112, 1761 (2014).
  * Tateishi  _et al._ [2014] K. Tateishi, M. Negoro, S. Nishida, A. Kagawa, Y. Morita, and M. Kitagawa, Room temperature hyperpolarization of nuclear spins in bulk, Proceedings of the National Academy of Sciences 111, 7527 (2014).
  * Niketic  _et al._ [2015] N. Niketic, B. Brandt, W. T. Wenckebach, J. Kohlbrecher, and P. Hautle, Polarization analysis in neutron small-angle scattering with a novel triplet dynamic nuclear polarization spin filter, Journal of Applied Crystallography 48, 1514 (2015).

Supplementary Information
Electron induced nanoscale nuclear spin relaxation probed by hyperpolarization
injection
William Beatrez,1 Arjun Pillai,1 Otto Janes,1Dieter Suter,2 and Ashok
Ajoy,1,3,∗

_1 Department of Chemistry, University of California, Berkeley, Berkeley, CA
94720, USA._

_2 Fakultät Physik, Technische Universität Dortmund, D-44221 Dortmund,
Germany._

_3 Chemical Science Division, Lawrence Berkeley National Laboratory,
University of California, Berkeley, Berkeley, CA 94720, USA._

Fig. S1: $T_{1}$ measurements.  Sample $T_{1}$ at low-field (36mT, blue) and
high-field (7T, red) at $\tau{=}60$s. Points are data and solid lines are
monoexponential fits. Dashed line denotes $1/e$, and intersections give a low-
field $T_{1}$ of 283$\pm$2s and a high field $T_{1}$ of 1520s$\pm$20s.  Fig.
S2: Polarization buildup with hyperpolarization period $\tau$. Panel shows the
unnormalized data corresponding to Fig. 3A of the main paper. The different
curves are the ${}^{13}\textrm{C}$ pulsed spin-lock decays obtained under
varying hyperpolarization time $\tau$ (see colorbar). This data is also shown
as a movie in Ref. [40].

##  I $T_{1}$ measurements at low and high field

Fig. S1 considers the measurement of the $T_{1}$ of the ${}^{13}\textrm{C}$
nuclear spins at the polarizing field ($B_{\textrm{pol}}{=}36$mT) and the
readout field ($B_{0}{=}7$T) respectively. For these measurements, the sample
is first hyperpolarized, and is subsequently carried to the field of interest
($B_{\textrm{pol}}$ or $B_{0}$) for a waiting period $t_{0}$. Finally, the
sample is transported to $B_{0}$ where spin-lock readout is carried out
similar to Fig. 2, and the integrated signal is plotted. The rapid field
cycling for these measurements is carried out via mechanical sample shuttling,
and the experimental strategy is similar to that in Ref. [27].

Notably, varying the period $t_{0}$ allows one to discern the profile of the
$T_{1}$ relaxation at both fields (see Fig. S1). We find that the relaxation
profiles closely follow monoexponential decay profiles (solid lines). From a
$1/e$ intercept (dashed line) in Fig. S1 we find the respective low and high
field relaxation times as $T_{1,\textrm{LF}}{=}283{\pm}$2s and
$T_{1,\textrm{HF}}{=}1520\pm 20$s respectively. We note that the
monoexponential decays here stand in contrast to the observed stretched
exponential decays that we find for the rotating frame lifetimes
$T_{2}^{\prime}$ in Fig. 2 and Fig. 3 of the main paper.

##  II Movies showing data in Fig. 2 and Fig. 3

As a complement to the graphs in the main paper, we present movies
corresponding to the data in Fig. 2A and Fig. 3A on Youtube (found at Refs.
[38, 39]). These movies show clearly the progressive slowing down of the decay
profiles upon increasing the hyperpolarization time $\tau$ in both
representations. The gray lines here show fitted stretched exponential lines
corresponding to the previous data for clarity, allowing a guide to the eye to
track the slowing relaxation dynamics with increasing $\tau$. Fig. S2 shows an
unnormalized plot of Fig. 3, from where the polarization buildup dynamics as a
function of $\tau$ can be extracted.


# Hierarchical Automatic Multilayer Power Plane Generation with Genetic
Optimization and Multilayer Perceptron

HAIGUANG LIAO1 VINAY PATIL1 XULIANG DONG1 DEVIKA SHANBHAG1
ELIAS FALLON2 TAYLOR HOGAN2 MIRKO SPASOJEVIC2 LEVENT BURAK KARA1
1\. Carnegie Mellon University
Pittsburgh, PA 15213
2\. Cadence Design Systems
San Jose, CA 95134 Address all correspondence to this author.

###### Abstract

We present an automatic multilayer power plane generation method to accelerate
the design of printed circuit boards (PCB). In PCB design, while automatic
solvers have been developed to predict important indicators such as the IR-
drop, power integrity, and signal integrity, the generation of the power plane
itself still largely relies on laborious manual methods. Our automatic power
plane generation approach is based on genetic optimization combined with a
multilayer perceptron and is able to automatically generate power planes
across a diverse set of problems with varying levels of difficulty. Our method
GOMLP consists of an outer loop genetic optimizer (GO) and an inner loop
multi-layer perceptron (MLP) that generate power planes automatically. The
critical elements of our approach include contour detection, feature
expansion, and a distance measure to enable island-minimizing complex power
plane generation. We compare our approach to a baseline solution based on A*.
The A* method consisting of a sequential island generation and merging process
which can produce less than ideal solutions. Our experimental results show
that on single layer power plane problems, our method outperforms A* in 71% of
the problems with varying levels of board layout difficulty. We further
describe H-GOMLP, which extends GOMLP to multilayer power plane problems using
hierarchical clustering and net similarities based on the Hausdorff distance.

##  1 Introduction

Figure 1: Single layer power plane generation on a PCB board. Middle: Pins
belonging to the same net are colored similarly. Right: Planes belonging to
the same net are colored similarly. Note that some nets are split into
multiple islands.

Following Moore’s Law [1], modern electronics designs are getting increasingly
more complex. Printed Circuit Boards (PCB) design, as one of the critical
steps in electronics design is also getting increasingly more complex and
time-consuming with the increasing size of boards, more complex design rule
constraints (DRCs), as well as signal integrity (SI) and power integrity (PI)
requirements. One of the most critical steps that hamper PCB design
development is the generation of a power plane [2], which aims to configure
the various metal planes’ layout on a PCB. Although there has been significant
development in automatic solvers for the critical objectives of PCB design
related to power plane generation such as IR-drop solver [3, 4, 5], power
integrity solver [6, 7], and signal integrity [8, 9, 10], the generation of
power plane itself still relies on the manual effort of electronics designers.
Thus, automatic power plane generation is needed to address an unmet need in
PCB design.

In this work, we present an evolutionary automatic power plane generation
method based on a multilayer perceptron (MLP) and genetic optimization (GO),
which we refer to as GOMLP. This method is based on a nested combination of a
multilayer perceptron and genetic optimization, with critical components
including contour detection, feature expansion and customized distance
metrics. Systematic experiments are conducted to study the contributions of
the critical components and also stress test the method on a diverse set of
problems with varying difficulties. The method is compared against a baseline
solution based on A* [11]. To the best of our knowledge, this work is the
first effort to solve power plane generation with fully automated methods. The
GOMLP automatic power plane generation method is then further extended to
solve multilayer power plane generation problem by integrating it with a
hierarchical method that solves the challenging combinatorial problem of layer
assignment. The hierarchical method applies hierarchical clustering with
customized distance metrics based on Hausdorff Distance and Earth Mover
Distance to generate dendrograms of nets. The generated dendrogram is applied
to determine layer assignment of nets in flexible and efficient ways. We refer
to this approach as Hierarchical GOMLP (H-GOMLP).

##  2 Background

###  2.1 Problem Formulation: Power Plane Generation

Power plane generation is a critical step PCB design, where the layout of
metal planes connecting critical components is configured on each layer of the
PCB board [12, 2]. As shown in Fig 1, a layer of the PCB board contains a set
of $m$ nets $\\{N_{1},N_{2},N_{3},...,N_{m}\\}$. Each net is composed of a
potentially large number of pins (the I/O ports of the components on the
board) $P_{N_{i}}=\\{p_{1},..,p_{q_{i}}\\}$, where $q_{i}$ is the total number
of pins for $N_{i}$. $p_{j}\in\mathbb{R}^{2}$ encodes the $x$ and $y$
coordinates of a pin. Net information and corresponding pin locations are
determined in a previous step of board design and component placement. As
such, all pin locations are fixed and are not subject to change. The objective
is to group together all the pins with the same label (net assignment) through
a metal plane, while excluding all other pins.

Power plane generation can be formulated as a mutually exclusive space
partitioning of the board $\Omega$ into $m$ distinct power planes
$\Omega_{1},..,\Omega_{m}$ such that:

$\Omega=\Omega_{1}\cup\Omega_{2}\cup...\cup\Omega_{m}~{}~{}and~{}~{}\Omega_{i}\cap\Omega_{j}=\emptyset$
(1)

Plane $\Omega_{i}$ corresponds to net $N_{i}$, which means the number of
planes is equivalent to the number of nets. Each plane $\Omega_{i}$ may
consists of several islands: $\Omega_{i}=I_{i}^{1}\cup I_{i}^{2}...\cup
I_{i}^{s_{i}}$, where $s_{i}$ is the total number of islands of plane (or net)
$i$.

Design Constraint: All pins of Net $N_{i}$ must be located inside $\Omega_{i}$
(possibly spread over multiple islands) and no pin from another net can appear
in $\Omega_{i}$:

$P_{N_{i}}=\\{p_{1},..,p_{q_{i}}\\}\in\Omega_{i}$ (2)

Design Objective: The power plane generation problem formulation described
above is the basic form. For each plane $\Omega_{i}$, it is preferable to
generate a singly connected plane consisting of a single island. Otherwise,
disconnected islands belonging to the same plane will have to be connected on
a separate layer of the PCB using extraneous connections and vias.
Mathematically, this can be stated as finding an optimized space partitioning
$\Omega^{*}$ which minimizes the total number of split islands:

$\Omega^{*}=argmin_{\Omega}\sum_{i=1}^{m}s_{i}$ (3)

Given a board, we denote a set of power planes that satisfy the design
constraint as feasible, and for each net, those that have a minimum number of
islands as desirable. Fig. 3 shows examples. For the example design, there are
two nets $\\{N_{1},N_{2}\\}$, where $N_{1}$ has two pins and $N_{2}$ has one
pin. The desirable designs (E and F) are with singly connected partitions that
are feasible. Designs C and D, while feasible, split one of the nets into two
separate islands. Compared to E and F, these designs are not desirable.

The objective of minimizing the number of islands for each net makes the
problem non-trivial since planes corresponding to different nets have to
compete for the shared design space while satisfying the design constraint.
Intuitively, single layer power plane generation is analogous to a 2D multi-
class image segmentation problem aimed at minimizing the number of islands for
all nets while ensuring that no pin is misclassified. While approaches such as
image segmentation [13, 14, 15] and space partitioning [16, 17] have been
extensively studied, these approaches do not put emphasis on split island
minimization. To the best of our knowledge the power plane generation is an
open problem not only in the scope of PCB design but also in the broader
context of space partitioning problems.

Multilayer Formulation: Most industry PCB designs consist of multiple layers
of boards and thus power plane generation solutions need to be extended to
account for multiple layers. Fig. 2 shows our workflow for multilayer cases.
Compared to the single layer case as shown in Fig. 1, the multilayer case
deals with one additional layer assignment step, where given a prescribed
number of layers $K$, each net in a design $N_{i}$ is assigned to a unique
layer $L_{i}$. After layer assignment, each layer $L_{i}$ consists of a set of
$z$ nets: $\\{N_{L_{i}1},N_{L_{i}2},...,N_{L_{i}z_{L_{i}}}\\}$. At this point,
the problem becomes a set of single layer problems where the power planes for
each layer are generated individually. The layer assignment step is
challenging owing to its combinatorial nature: Enumerating all possible
choices of layer assignment with $N$ nets and $K$ layers requires solving
single layer power plane problems $N^{K+1}$ times, which is intractable with
limited computation resources. This motivates us to develop an efficient and
flexible hierarchical method to solve multilayer power plane generation
problems.

Figure 2: Multilayer power plane generation on a PCB board.

###  2.2 Scope

Our approach seeks to find a feasible solution that minimizes the number of
islands. In this work, the priority is on identifying power planes that match
this desire. However, we do not optimize the outer contour geometries of the
planes and reserve this consideration for future work. As such, in the
presented work, designs E and F in Fig. 3 are identical in quality.

In industrial PCB design, there exists other constraints and objectives for
space partitioning. These include considerations around IR-drop, power
integrity, and signal integrity of the partitioned board. For instance, it is
common to impose an upper bound on the admissible IR-drop for each source-
drain pin pair within a plane, which favors planes with the least number of
thin regions. For each power plane candidate, the evaluation of the IR drop
involves a physical simulation that can be prohibitive. In this work, we thus
exclude such considerations in favor of establishing the foundations for
island minimizing partitioning. Future work will have to account for such
criteria.

In this work, GOMLP and H-GOMLP are used to solve different power plane
problems. GOMLP only considers a single layer of PCB, aiming to identify an
optimized solution to the nets and pins already assigned to a specific layer
(i.e., net and pin assignment is assumed to have been completed). By contrast,
H-GOMLP, which extends GOMLP, focuses on solving multilayer power plane
generation. In general, while GOMLP can be viewed as a component of H-GOMLP,
it can be utilized as a standalone, single layer PCB generation tool.

Figure 3: Plane feasibility vs. desirability.

##  3 Methods

Figure 4: Proposed GOMLP method to solve single layer power plane generation
problems.

Since GOMLP is an independent module but also a component of the more complex
H-GOMLP, we will first introduce GOMLP as an standalone method and then
introduce H-GOMLP.

In GOMLP, we parameterize the problem using a set of handles defined for each
plane $\Omega_{i}$ as $H_{N_{i}}=\\{h_{1},..,h_{k_{i}}\\}$, where
$h_{j}\in\mathbb{R}^{2}$ encodes the $x$ and $y$ coordinates of the handle on
the board. These handles serve as labelled movable points on the board such
that an optimized placement of the entire set of handles $H_{N_{i}},i:1..m$
produce a partitioning of the board where all original (and static) pin sets
$P_{N_{i}},i:1..m$ are correctly grouped, and the number of power plane
islands is minimized.

GOMLP is summarized in Fig. 4. It consists of an outer loop genetic optimizer
(GO) that identifies the best handle placements. To evaluate a candidate
solution, we use an inner loop multi-layer perceptron (MLP) that produces a 2D
segmentation of the board using $H_{N_{i}}$ and $P_{N_{i}},i:1..m$ as the
labelled data to be classified. The MLP uses an expanded feature space and, is
allowed to overfit to the data so as to ensure all pins are correctly
classified (i.e., satisfying the design constraint). A separate module then
identifies the total number of islands together with the pairwise distances
between the islands of the same net. These quantities are appended to GO’s
fitness function to drive the optimal handle placement.

###  3.1 Genetic optimization for handle placement

Genetic Optimization (GO) has found extensive use in problems where the
gradients of the objective function and constraint are not readily available
[18, 19, 20]. In this work, we use GO to optimize the handles’ locations,
which then are used as input to the MLP to produce the space partitions . The
primary motivation for applying GO instead of other gradient-based optimizers
is that the objective function in the power plane generation problem is not
readily differentiable.

In our GO, a population consists of $M$ chromosomes, where each chromosome
encodes the handle coordinates as a vector of genes. Thus, each chromosome is
a candidate solution for $\Omega$. Each gene has two scalars $x,y$ encoding
the 2D coordinates of a handle. For each $\Omega_{i},i:1..m$, we instantiate
$k$ handles, where k is calculated as follows:

$k=\frac{2\cdot\sum_{i=1}^{m}q_{i}}{m}$ (4)

where $q_{i}$ is the total number of pins for net $i$. Each chromosome is thus
a $2\cdot k\cdot m\times 1$ real vector, with the handles belonging to each
net occupying a $2\cdot k\times 1$ block. Handles are randomly initialized in
the first generation of GO. We use a fitness-proportional parent selection,
uniform mutation and probabilistic random swap for crossover operators [21,
22]. Fig. 4-top shows the collection of pins (letters) and handles (circles)
for a 5-net design problem.

###  3.2 Power plane generation using MLP

Given a chromosome, we use an MLP to construct the partitioned power planes.
The MLP takes as input the labelled training set $P_{N_{i}}\cup
H_{N_{i}},i:1..m$ and creates a region segmentation of the board using the
multi-class cross entropy loss. Fig. 4-middle shows the current space
partitioning result. For each item in the training set, we expand its $x,y$
coordinates into a $15\times 1$ vector: $\\{x,y,x\cdot y,\sin(2\pi\cdot
x),\cos(2\pi\cdot x),x^{2},\sin(2\pi\cdot y),\cos(2\pi\cdot
y),y^{2},\sin(3\pi\cdot x),\cos(3\pi\cdot x),x^{3},\sin(3\pi\cdot
y),\cos(3\pi\cdot y),y^{3}\\}$ to enable complex boundaries. The MLP is
trained with these features as input and the net classification as output is
then used to generate space partition $\Omega$ by predicting the plane
assignment for every coordinate from the meshgrid on the PCB board layer.

The MLP consists of three hidden layers, with 50 neurons in each layer and the
$tanh$ activation function. Adam optimizer with a learning rate of 0.002 is
used. It is desirable to overfit to the training data, so as to favor all pins
to satisfy the feasibility constraint of Eqn. 2. Nonetheless, the resulting
partitioning may result in disconnected islands, and may not satisfy Eqn. 2.

###  3.3 Fitness calculation

To evaluate each solution’s fitness, we determine the total number of islands
$s_{i}$ of each plane $\Omega_{i}$ generated by the MLP using a contour
detection algorithm [23].

$F^{island}=-\sum_{i=1}^{m}s_{i}$ (5)

Additionally, we append the fitness function with a measure to distinguish the
islands of a plane $\Omega_{i}$ that are far apart from those that are
proximate. This helps GO acquire guidance in handle placement for the
subsequent generation (Fig. 4-bottom). For $\Omega_{i}$, we define:

$F_{i}^{d_{min}}=-\sum_{a,b}^{s_{i}}||I_{i}^{a}-I_{i}^{b}||_{min}$ (6)

where $I_{i}^{a},I_{i}^{b}$ denote islands $a$ and $b$ belonging to
$\Omega_{i}$, and $||.||$ is the Euclidean norm. $F_{i}^{d_{min}}$ captures
the sum of minimum distances between the pairs of the disjoint islands of
$\Omega_{i}$.

While min. distances are useful in signaling nearly connected islands, they
typically favor the generation of elongated thin connections. To additionally
drive the islands toward one another, we define a centroid distance as a part
of the fitness function:

$F_{i}^{d_{cent}}=-\sum_{a,b}^{s_{i}}||I_{i}^{a}-I_{i}^{b}||_{cent}$ (7)

Fig. 5 illustrates the above two islands distance calculations. For a
segmentation solution produced by the MLP, the final fitness maximized by the
GO is:

$F=F^{island}+\sum_{i}^{m}F_{i}^{d_{min}}+\sum_{i}^{m}F_{i}^{d_{cent}}$ (8)

Figure 5: Contour detection and distance metric calculations of space
partitions. Figure 6: Details of GO structures and flow

Input : Netlist of designs: $\\{N_{1},N_{2},...,N_{m}\\}$

Output : Space partitions: $\\{\Omega_{1},\Omega_{2},...,\Omega_{m}\\}$

Initialize handles $\\{h_{1},h_{2},...,h_{k_{i}}\\}$ for each net $N_{i}$ as
the first population for GO;

for  _generation=1,…, $N$_ do

MLP fitting all fixed points $\\{P_{1},P_{2},...,P_{q_{i}}\\}$ and handles
$\\{H_{1},H_{2},...,H_{k_{i}}\\}$ for each net $N_{i}$ with augmented
features;

Contour detection and distance metric calculations for all planes
$\\{\Omega_{1},\Omega_{2},...,\Omega_{m}\\}$ ;

Calculating fitness scores for all chromosomes based on MLP fitting results;

Generating new population with crossover and mutation based on selected
elites;

end for

Algorithm 1 GOMLP power plane generation

###  3.4 H-GOMLP

H-GOMLP extends GOMLP to multilayer power plane generation problems. Fig. 7
shows the flowchart of our proposed H-GOMLP method, which consists of 3 major
steps:

  * •

1\. Calculating pairwise distances between all net pairs

  * •

2\. Running hierarchical clustering for all nets to facilitate layer
assignment

  * •

3\. Running GOMLP on each layer

In the first step, given all nets $\\{N_{1},N_{2},...,N_{m}\\}$ in a
multilayer PCB board to be designed. We calculate the Hausdorff Distance (HD)
[24] and Earth Mover Distance (EMD) [25] between all pairs of nets. For a pair
of nets $\\{N_{i},N_{j}\\}$, the distance is calculated based on the
coordinates of their corresponding pins. For HD, the distance is calculated
by:

$\begin{split}d_{HD}(N_{1},N_{2})=max\\{max_{x_{i}\in N_{1}}min_{y_{j}\in
N_{2}}||x_{i}-y_{j}||,\\\ max_{y_{j}\in N_{2}}min_{x_{i}\in
N_{1}}||x_{i}-y_{j}||\\}\end{split}$ (9)

For EMD, the distance is calculated by:

$\begin{split}d_{EMD}(N_{1},N_{2})=min_{\phi:N_{1}\rightarrow
N_{2}}\sum_{(x,y)\in N_{1}}||(x,y)-\phi((x,y))||_{2}\end{split}$ (10)

where $\phi:N_{1}\rightarrow N_{2}$ is a bijection. The motivation for using
the HD and EMD to calculate the distance between each pair of nets is to
approximate the proximity and overlapping of the nets’ layouts. Smaller HD and
EMD distances indicate that the pins’ layouts of two nets are very close to
each oother or significantly overlap. In this way, we can make use of this
information in the layer assignment step as an effective metric for the
spacing of nets.

After the distance calculation, a hierarchical clustering for all the nets is
done based on the pairwise distance of all nets. In the hierarchical
clustering, instead of using the typical distance metrics like Euclidean
distance, we use the inverse of the HD and EMD between pairs of nets
calculated in previous step. The motivation of using inverse of HD and EMD is
to cluster nets that are less overlapping on to the same layer in order to use
each layer most efficiently. In such case, the relative difficulty of
finishing up those single layers with GOMLP would be minimum. In another
respect, the overall usage of each layer space is more efficient.

The hierarchical clustering will generate a dendrogram of nets as shown in
Fig. 7. Based on this dendrogram, we can generate clusters of nets with
specified budget of cluster numbers. Nets belonging to the same cluster will
then be assigned to a unique layer of the PCB board. As shown in Fig. 7, based
on the dendrogram, we can get a 2-layer case at a higher level on the
dendrogram and get another 3-layer case at a slightly lower level. This is a
unique advantage of our approach, since it provides flexibility in the number
of layers needed for a design: for a given multilayer PCB board, the number of
layers needed to finish power plane generation is usually not given upfront.
With our approach, the designer can start at a very high level on the
dendrogram, which gives a small layer number solution. If the small layer
number does not lead to feasible or desirable designs, the designer can move
down on the dendrogram to get another design with incrementally more number of
layers, and the process can be repeated until desirable or feasible power
planes are generated.

Finally, after the layer assignment, the H-GOMLP is finished by running GOMLP
for all the layers separately and the details of GOMLP have been described
above. This step and the previous layer assignment step can be executed
iteratively if necessary. For instance, after the first run of layer
assignment with total number of layers determined to be $K$, GOMLP is applied
to solve those $K$ layers independently. If it turns out that some layers
cannot achieve desirable or feasible designs with GOMLP, we can go back to the
layer assignment step and generate a new assignment with $K+1$ layers,
followed by running GOMLP for these $K+1$ layers. One great advantage for our
approach is that GOMLP execution for different layers can run in parallel as
well as the GO module inside GOMLP loop, which is already realized in our
implementation. This enables our method’s scalability to large scale industry
PCB board designs. Detailed description of H-GOMLP is available in Algo. 2.

Figure 7: Proposed H-GOMLP method to solve multilayer power plane generation
problems.

Input : Netlist of designs: $\\{N_{1},N_{2},...,N_{m}\\}$

Output : Space partitions for $K$ layers:
$\\{\Omega_{i1},\Omega_{i2},...,\Omega_{in}\\}$, where $i=1,...,K$

Calculate the HD and EMD for all pairs of nets $\\{N_{i},N_{j}\\}$;

Hierarchical clustering and generate dendrogram for all nets based on the
inverse of HD/EMD;

Layer assignment for all nets with specific number of layers $K$ based on the
dendrogram;

for  _layer=1,…, $K$_ do

Run GOMLP for nets assigned on the layer to generate space partitions:
$\\{\Omega_{i1},\Omega_{i2},...,\Omega_{in}\\}$

end for

Algorithm 2 H-GOMLP power plane generation

##  4 Method: A* Baseline Method

A* is one of the most frequently used search algorithms in graph traversal and
path finding. In this work, A* search is used as the backbone of our baseline
power plane generation algorithm, which we refer to as A*. It has been proven
that with admissible heuristics, the A* can find the shortest path between two
vertices in a graph [26]. Fig. 8 shows the flowchart of the A* baseline. The
reason for selecting the A* as our baseline is that it has been the most
frequently used algorithm practice to solve routing problem [27, 28, 29]. In
the baseline method, the solution consists of multiple steps and one of the
key steps (Island Connection) can be formulated as a routing problem.

The input of the power plane generation problem is clusters of pins belonging
to each net $\\{N_{1},N_{2},...,N_{m}\\}$. For A*, the first step is the
Island Generation. At this step, for each net $N_{i}$, a minimum spanning tree
($MST_{N_{i}}$) among all the pins $P_{N_{i}}$ belonging to net $N_{i}$ is
generated with using Kruskal’s algorithm [30]. After the MST for all nets
$\\{MST_{N_{1}},MST_{N_{2}},...,MST_{N_{m}}\\}$ are generated, edge pruning is
applied to all the MSTs where any edges of the MSTs that intercept edges of
other nets are pruned. After pruning, there is a set of disconnected trees
(pins connected with edges) for each net, which is represented by an island
that includes all the pins belonging to the tree while excluding all other
pins not belonging to the tree. For instance, if island generation produces
$s_{i}$ disconnected trees $\\{T_{1},T_{2},...,T_{s_{i}}\\}$ for Net $N_{i}$,
then the disconnected trees are represented by a corresponding set of islands
$\\{I^{1}_{i},I^{2}_{i},...,I^{s_{i}}_{i}\\}$, where the area of $I_{j}$
include all the pins of tree $T_{j}$ but not include pins belonging to other
trees of net $N_{i}$ or other nets. Based on the Island Generation step, each
net is represented as a set of islands. The islands for all the nets are fed
into the downstream Island Connection step, where A* algorithm is applied to
route all the disconnected islands for each net to generate a connected tree
for each net. Specifically, given an order of the nets, for each net $N_{i}$
with disconnected island set $\\{I^{1}_{i},I^{2}_{i},...,I^{s_{i}}_{i}\\}$, A*
works on connecting the disconnected island set. For the disconnected island
set with $s$ islands, $(s-1)$ pairs of islands are generated, and then A* will
sequentially connect all the pairs of islands by searching through a graph
consisting of upsampling nodes distributed on the board. The Island Connection
is finished when all the island pairs belonging to all the nets have been
attempted to be connected with A*.

Since all the island pairs of one design have to share the same space when
connected with A*, it is not guaranteed that all of them can be connected. If
the Island Connection step is successful, then a set of trees
$\\{T_{N_{1}},T_{N_{2}},...,T_{Nm}\\}$ will be generated where $T_{N_{i}}$
connected all the pins belonging to net $N_{i}$. If some islands are not fully
connected, some nets will end up having more than one tree. Based on the
generated tree set from the Island Connection step, the final power plane is
generated by inflating the tree edges of each tree. Specifically, densely
located nodes are sampled on all edges of the trees so that each tree
$T_{N_{i}}$ can be represented as a set of nodes. Then the space partition is
generated by running K-nearest Neighbors (KNN) [31] algorithms based on the
nodes from all trees with the number of neighbor $k$ set as 1.

Figure 8: Flowchart of A* Power Plane Generation Method.

The A* baseline approach provides a practical solution flow to generate a
power plane. However, a downside of A* is its greedy manner [32, 33]: it
focuses on finding out the optimal solution for only one search problem. In
power plane generation, there are cases where multiple search problems that
need to be solved on the same graph sharing limited graph space. Under such
case, A* would struggle to consider all the problems simultaneously. Thus,
there are at least the following combinatorial optimization challenges to be
addressed to make it an applicable approach:

  * •

Sequence of nets: in the Island Connection step, the A* has to be applied to
the nets sequentially. Thus, an optimal sequence of nets needs to be
determined. Otherwise, the method will be unfavorable to the nets that are
routed later because the nets routed earlier have occupied part of the search
graph.

  * •

Sequence of island pairs: similar to the sequence of nets challenge, within
each net, the island pairs also need to be sequenced before A* is applied. The
sequence is a non-trivial combinatorial optimization problem.

  * •

Design of upsampling nodes: to run the A* algorithm efficiently, the search
graph size need to be reasonably small, and thus designing the distribution of
upsampling nodes is another combinatorial optimization challenge that needs to
be solved.

Input : Netlist of designs: $\\{N_{1},N_{2},...,N_{m}\\}$

Output : Space partitions: $\\{\Omega_{1},\Omega_{2},...,\Omega_{m}\\}$

Tree generatoin $\\{MST_{N_{1}},MST_{N_{2}},...,MST_{N_{m}}\\}$ with MST for
each net;

Island generation by tree pruning for each net $N_{i}$:
$\\{T_{1},T_{2},...,T_{s_{i}}\\}$;

for net for  _i=1,…, $m$_ do

Running A* to connect the isolated islands
$\\{I_{i}^{1},I_{i}^{2},..,I_{i}^{s_{i}}\\}$ for net $N_{i}$, forming complete
tree edges for each net

end for

Generating final space partitions by running KNN on the upsampled nodes from
edges of different net

Algorithm 3 A* power plane generation

##  5 Evaluation

To systematically study the performance of our approach, extensive experiments
are performed, which can be divided into two parts. In part one, GOMLP is
applied to solve problems with and without the key components to show the
effect of the key features of our method, which include the following:

Effect of GO on MLP: Given the same design problem, power planes are generated
with standalone MLP and MLP with GO to investigate the effect of GO upon MLP.

Effect of feature expansion on GOMLP : Given the same design problem, GOMLP is
applied with and without expanded features to investigate the effect of
feature expansion.

Effect of distance metric on the GOMLP: Given the same design problem, GOMLP
is applied with and without the distance metric to investigate its effects on
performance of the GOMLP.

In the second part of the experiments, we apply GOMLP as well as the baseline
A* method to a set of power plane generation problems with different
configurations and difficulties to systematically study the performance of our
method.

###  5.1 Dataset

To generate the problems, given a layer of design on a PCB board, we vary the
number of nets from 6 to 8, with more nets corresponding to more challenging
problems. For each net number, combinations of nets are selected to generate
distinct problems. With the variations, there are altogether 129 different
problems with varying difficulties. The GOMLP as well as the baseline A* are
applied to each of the problems with maximum run time for solving each problem
set as 30 minutes for a fair comparison.

Cross entropy is used as the loss function. For the GO, the maximum number of
generations is set as 20, population size is 30 and elite size is 10. An
industry real PCB board design Beaglebone is used, and different design
problems in the experiment come from different design problems on different
layers. Since the objective in this work is to minimizes the total number of
split islands, in evaluating the quality of generated power planes, the
primary metric used is Extra Islands (EI), which is defined as follows:

$EI=\sum_{i=1}^{m}s_{i}-m$ (11)

For evaluating H-GOMLP method on multilayer PCB board designs, we generate
more challenging problems by increasing the number of nets to 20, which is
also based on industry real design Beaglebone. We set the parameters of GOMLP
inside H-GOMLP the same as the stand-alone GOMLP. In order to compare
effectiveness of different settings of H-GOMLP, the primary metric we use is
the minimum number of layers we need to finish a multilayer power plane
generation problem where all layers have 0 EIs. We call this metric: Minimum
Close Design Layer (MCDL). It is worth mentioning that the MCDL would change
based on different parameter settings of GOMLP or difficulty level of the
designs and thus the conclusion based on MCDL only holds for specific problems
and settings.

##  6 Results and Discussions

###  6.1 Comparison between MLP and GOMLP

Figure 9: Power plane generation on 5 nets PCB layer design case using MLP vs
GOMLP. MLP: 6 partitions with 2 disjoint partitions belonging to net E ,
evolutionary MLP : 5 partitions with no disjoint partitions

Fig. 9 , shows the results of MLP and GOMLP on a 5 net case. On the left, the
partitions are generated using only MLP. Although the MLP can generate a power
plane design that correctly classifies the nets. It fails to converge to the
most desired power plane partitioning where the partition for each net is
connected. As a result, the partition corresponding to net E is divided into
two islands. On the right side, the partition is generated using GOMLP. It was
executed for just 2 generations and generated the desired result which
satisfied constraints. Also, more importantly, the space partitions for each
net are all single connected planes. The result demonstrates that the power
plane generation is a nontrivial space partition problem that can not be
solved directly with a classifier. The challenge for a single classifier to
generate a feasible power plane is due to the difficulty of encoding the
design objective of overall island minimization as a differentiable function
into the loss function of MLP. While Euler number [34] allows the
determination of the number of disconnected islands, connecting that loss to
the handle parameters is not readily achievable.

In our GOMLP, the combination of MLP and GO solves this challenge by
optimizing the handle locations to control the shape of the MLP generated
space partitions. Since handle optimization seeks to minimize the number of
islands, GO will work on connecting isolated islands of the MLP-generated
partitions implicitly.

###  6.2 Comparison between GOMLP with and without distance based cost

Fig. 10 shows the experiment results of the GOMLP with and without the
distance metric described in the Method part. A 6-net PCB board layer design
is given and the GOMLP is applied twice to solve the problem: the first time
with the distance metric and the second time without the distance metric. The
distance metric describes the proximity of isolated islands for each net. On
the left of Fig. 10, the generated power plane is based on GOMLP without
including distance metric, which means only the number of isolated islands are
used to describe the layout of the space partition. On the power plane, the
space partition of net E is divided into two isolated islands on the left and
right sides of the plane. While on the right of Fig. 10, the generated power
plane is based on GOMLP with the distance metric included. In this case, the
space partition corresponding to net E is connected as a single U-shaped plane
during the iterative optimization procedure of the GOMLP.

The improved performance with the introduction of the distance metrics
compared to cases where the only number of isolated islands are used in cost
function can be explained as follows: in the optimization process of the
GOMLP, the space partition is controlled by the layout of the handles. The
handles are iteratively optimized to generate higher fitness scores through
the GO operation. This mechanism means that the GO is working on optimizing
the planes implicitly through the handles’ movement. As a result, in cases
where the handles for one net are distributed as distant clusters, for
instance, the case in Fig. 10, it would be challenging for the GOMLP to
provide the handles with a good moving direction without the introduction of
distance metric. In other words, the number of isolated islands is not enough
to provide enough guidance for driving the handles to gradually move and form
connected space partitions for each net. The distance metric provides
something similar to gradient information to guide the handles gradually move
in order to connect isolated islands for each net.

Figure 10: Power plane generation on 6 nets PCB layer design using GOMLP with
and without distance metric.

###  6.3 Comparison of GOMLP with and without feature expansion

Figure 11: Power plane generation on 6 nets PCB layer design using GOMLP with
and without feature expansion. GOMLP: 7 partitions with 1 disjoint partitions
belonging to net D, with feature expansion: 6 partitions with no disjoint
partitions.

Fig. 11 shows the results GOMLP with and without feature expansion on a 6 net
PCB layer design. On the left, the partitions are generated with GOMLP without
the use of feature expansion. Without the expanded features, the model is
unable to generate complex space partitions that are feasible and desirable.
In the generated power plane, $\Omega_{D}$ is separated into two isolated
islands, which makes it undesirable. Whereas on the right of Fig. 11, the
partitions are generated by GOMLP with feature expansion. With the additional
higher order features as input, the MLP can learn more complex space
partitions and curved space partitions which leads to desirable power plane
generation. In the generated power plane, the partition for $\Omega_{D}$ is a
singly connected plane with a thin channel connecting the two isolated regions
of net D. The comparative results demonstrate the effect of feature expansion
in our method. The introduction of high dimensional features makes it easier
for the MLP to generate complex space partitions such as zig-zag-shaped
planes.

###  6.4 Scalability of GOMLP

As described in the Experiment part, the second part of the experiments
includes the applications of both the GOMLP and A* method to problems with
different configurations and difficulties to systematically assess the
performance and scalability of the GOMLP and A*. For the 129 problems
experimented, the number of nets ranged from 6 to 8, and problems with more
tend nets problems tend to be more difficult. Given a solution, EI is used as
the primary metric. Only when EI equals 0, the generated power plane is
desirable.

Fig. 12 shows the results of GOMLP and A* on the set of 129 problems. The plot
shows the number of extra islands for solutions of GOMLP and A* in the
ascending order of evolutionary MLP results, as well as their gaps. The plot
shows that GOMLP is significantly better in generating high quality power
planes with fewer extra islands. Among all the 129 experiments, GOMLP
outperforms or has the number of extra islands as A* in 71% of cases. For
problems where A* and GOMLP perform similarly, both methods tend to identify
the optimal power planes ($EI=0$). A paired t-test of GOMLP and A* on our 129
board dataset indicates that GOMLP is statistically better than A* (p<0.001).

A* outperforms GOMLP in only 4 out of the 129 experiments. Detailed studies
are conducted on those rare cases where A* outperforms GOMLP. The primary
reason that GOMLP did not perform as well as A* in these cases is due to the
presence of very near-by pins that belong to different nets, challenging the
MLP to correctly classify these pins. However, these cases will not likely be
encountered in realistic PCB board designs given the multilayer flexibility of
the problem: designers would move some of the very close pins belonging to
different layers of the PCB board to avoid such corner cases. To sum up, GOMLP
has good scalability to non-trivial power plane generation designs and
performs significantly better than A* baseline.

Figure 12: GOMLP vs A* for the 129 different board problems. Vertical axis
shows the extra islands. GOMLP is better than A* in 71% of the cases. Figure
13: Resulting power planes for a 5-net design problem Figure 14: Evolution of
the handles across the generations. Colored circles: handle locations, letters
A-E: pins belonging to different nets. Elite$\\#4$ in generation$\\#3$
converges to the desired result. Figure 15: GOMLP (top) and A* (bottom)
generated power planes for problems with 6,7 and 8 nets.

Finally, to give a more intuitive demonstration of how the GOMLP works, Fig.
13 and Fig. 14 show the evolution of space partitions and the handles across
generations. In Fig. 13, the pins (annotated with alphabets) are correctly
classified in the very first few generations, however, the space partition for
net E is divided into two separated islands in the left and right side of the
board. As the evolutionary MLP keeps optimizing the space partitions, the two
islands for net E start to approach and then finally merged into a desired one
single plane. The handles are the implicit driving force for the evolution, as
shown in Fig. 14, where the handles, pins, and the planes are plotted
together. The handles locations changed across generations to implicitly
change the space partitions generated by the MLP. It is worth mentioning that
the moving patterns of the handles are not as intuitive as we expected. For
instance, during the merge of two islands for net E, the handles for net E on
two different islands are not directly moving towards the other islands. This
is because the relationship between handles and the space partitions is highly
nonlinear and it is part of our future work to further study the evolutionary
patterns of the handles. For design C,D,E and F in Fig. 14, they are all
converged solutions for the same problem, showing another advantage for GOMLP:
generating multiple solutions for the same problem. Fig. 15 shows the
generated power planes by GOMLP and A* on the same set of the problems with 6,
7, and 8 number of nets. For the three problems, GOMLP is able to generate
desirable power planes, while A* method generates power planes that have more
EI numbers.

###  6.5 H-GOMLP: effectiveness of layer assignment

To validate the effectiveness of our chosen metrics (EMD and HD) in the
hierarchical clustering, we studied the structure of dendrograms generated
during hierarchical clustering of the H-GOMLP, especially the sequence by
which different pairs of nets are merged on the dendrograms. Fig. 16 shows
dendrograms based on hierarchical clustering of 20 nets in a multilayer power
plane generation problem based on the inverse of HD and the inverse of EMD. In
the dendrogram based on the inverse of HD, net pair A and net pair B are
highlighted and their pins are plotted as scatter plots respectively. Net pair
A is the pair that is merged first on the dendrogram and net pair B is the
last pair merged. By comparing their pins’ scatter plots, it can be seen that
pins of net pair A are further away to than those of net pair B, which
corresponds to the fact that the HD of net pair A is larger than that of net
pair B. In the dendrogram based on the inverse of EMD, the net pair C
(highlighted) is the first pair merged on the dendrogram while net pair D
(highlighted) is the last pair merged. Their respective scatter plots of pins
also show similar trend as those of HD case: pins for net pair C are further
away and less overlapping than pins for net pair D. The results show the
effectiveness of our chosen metrics in prioritizing the merge of nets with
less overlapping area or are further away. Clusters of nets (I,II,III,IV) that
are merged at the higher level on the dendrograms are also plotted.

Figure 16: Dendrograms from hierarchical clustering of 20 nets in a multilayer
power plane generation problem based on 1/HD and 1/EMD. Pair A and Pair C are
net pairs that are the first pairs merged; pair B and pair D are net pairs
that are the last pairs merged; clusters I,II,III and IV are the clusters of
nets that are merged at the highest level on the dendrograms.

###  6.6 H-GOMLP: results on solving non-trivial multilayer problems

Fig. 17 shows the results of H-GOMLP on solving a 20-nets multilayer power
plane generation problem with HD as the distance metric. The dendrogram is
shown on the left and it is marked with the level corresponding to 5 layers.
The corresponding power planes for all 5 layers in the first and last
generation of GOMLP are shown on the right. With 5 layers, H-GOMLP is able to
generate desirable power planes for all 5 layers. We also ran the 4-layers
case, and one of the layer ended up having number of EIs larger than 0, which
makes the design not desirable. So the MCDL for this problem solved with
H-GOMLP (with HD) is 5. For comparison, we also show the results for the same
20-net problem by H-GOMLP based on EMD as distance metric in Fig. 18. With the
same number of layers, the EMD-based H-GOMLP can achieve desirable designs
(EIs=0) for 4 out of the 5 layers while leaving one layer having larger than 0
number of EIs. This indicates that the MCDL for this problem and H-GOMLP (with
EMD) is larger than 5. Also, the EMD-based H-GOMLP tends to generate less
balanced layer assignment by assigning only one net to Layer 5. The smaller
MCDL and more balanced layer assignment of HD-based H-GOMLP both suggest that
it is better than the EMD-based H-GOMLP for this problem.

Figure 17: H-GOMLP generated dendrogram and multilayer power planes for a
20-nets multilayer board design based on inverse of HD. Figure 18: H-GOMLP
generated dendrogram and multilayer power planes for a 20-nets multilayer
board design based on inverse of EMD.

###  6.7 Discussions

In comparing GOMLP and A*, one obvious advantage that contributes to the
better performance of GOMLP is that it can avoid the combinatorial challenges
that A* faces: net sequencing, sequence of islands pairs as well as the design
of upsampling nodes. It has been shown [35] that there is no heuristic for net
ordering that is universally better for all possible problems. This means that
methods like A* need to solve the NP-hard combinatorial problem. Even if
optimized net sequences or upsamplings are given, the sequential nature of
such a method means that the subproblems solved first are prioritized over
those solved later, which is not desirable in partitioning problems of this
sort [36]. GOMLP provides a good alternative to avoid such combinatorial
challenges. Another advantage of the GOMLP over A* is the flexibility to add
complex objectives into the method. The A* method is based on generating
skeletons of each net using search, which solely focuses on finding the
shortest paths, this makes adding more complex objectives such as IR-drop to
A* challenging. On the other hand, GOMLP optimizes over the fitness score in
GO, which means the objectives do not have to be explicit and differentiable.
This will make sure that the GOMLP can be easily extended to more complex
objectives.

H-GOMLP is able to extend GOMLP to solve more challenging multilayer power
plane generation problems. One major advantage of H-GOMLP method comes from
its merge-based mechanism. In multilayer PCB board designs, one major
challenge is to determine the number of layers needed upfront. In most cases,
designers have to rely on tedious trial-and-error approaches which require
repeating of the time-consuming multilayer power plane generations with
different number of layers, due to lack of proper method to estimate the layer
numbers needed. The merge-based nature of H-GOMLP makes the multilayer power
plane generatioin a more natural bottom-up procedure, where the number of
layers does not need to be determined upfront. Instead, designers can easily
get layer assignment choices with incrementally more layers simply by getting
clusters from different levels of the dendrogram. It is worth mentioning that
for this problem, we don’t give the exact optimal power plane generation
solution for this 20-nets multilayer power plane generation because of the
impractically high computation required. For instance, given a layer number of
5, we need to run GOMLP for each layers of $20^{5}$ (3.2 million) combinations
of layer assignment cases. The optimality gap of H-GOMLP and further
improvement of H-GOMLP will be the focus on our future work.

##  7 Conclusion

In this work, we present a GOMLP automatic power plane generation method based
on the combination of MLP and GO. The method is compared against a baseline A*
solution. The GOMLP is based on a combination of multilayer perceptron and
genetic optimization, with critical components including contour detection,
feature expansion as well as the customized distance metric. The GOMLP
demonstrates the ability to automatically generate power planes across a
diverse set of non-trivial problems with different levels of difficulty.
Comparative study also shows the GOMLP is significantly better (in 71% cases)
than A* and those critical components including feature expansion, distance
metric as well as the combination of multilayer perceptron and genetic
optimization contribute to the competitive capabilities of the GOMLP in
automatically generating power planes. We further extend the GOMLP to H-GOMLP
which can solve multilayer power plane generation problems by combining GOMLP
with hierarchical clustering and customized distance metrics. Future works
include adding more complex objectives such as IR-drop to the GOMLP and
further improve the performance and scalability of H-GOMLP.

## Acknowledgment

This work is partially funded by the DARPA IDEA program (HR0011-18-3-0010).

## References

  * [1] Schaller, R. R., 1997.  “Moore’s law: past, present and future”.  IEEE spectrum,  34(6), pp. 52–59.
  * [2] Smith, L. D., Anderson, R., and Roy, T., 2001.  “Power plane spice models and simulated performance for materials and geometries”.  IEEE Transactions on Advanced Packaging,  24(3), pp. 277–287.
  * [3] Zhong, Y., and Wong, M. D., 2005.  “Fast algorithms for ir drop analysis in large power grid”.  In ICCAD-2005. IEEE/ACM International Conference on Computer-Aided Design, 2005., IEEE, pp. 351–357.
  * [4] Nithin, S., Shanmugam, G., and Chandrasekar, S., 2010.  “Dynamic voltage (ir) drop analysis and design closure: Issues and challenges”.  In 2010 11th International Symposium on Quality Electronic Design (ISQED), IEEE, pp. 611–617.
  * [5] Zhang, M.-S., and Tan, H. Z., 2015.  “Ir-drop modeling and reduction for high-performance printed circuit boards”.  IEEE Electromagnetic Compatibility Magazine,  4(4), pp. 90–101.
  * [6] Wu, T.-L., Chuang, H.-H., and Wang, T.-K., 2010.  “Overview of power integrity solutions on package and pcb: Decoupling and ebg isolation”.  IEEE Transactions on electromagnetic compatibility,  52(2), pp. 346–356.
  * [7] Fizesan, R., and Pitica, D., 2010.  “Simulation for power integrity to design a pcb for an optimum cost”.  In 2010 IEEE 16th International Symposium for Design and Technology in Electronic Packaging (SIITME), IEEE, pp. 141–146.
  * [8] Kim, J., and Li, E., 2010.  “Special issue on pcb level signal integrity, power integrity, and emc”.  IEEE Transactions on Electromagnetic Compatibility,  52(2), pp. 246–247.
  * [9] Chen, X.-P., 2010.  “Analysis and application for integrity of pcb signal”.  In 2010 2nd IEEE International Conference on Information and Financial Engineering, IEEE, pp. 328–331.
  * [10] Eudes, T., Ravelo, B., and Louis, A., 2011.  “Experimental validations of a simple pcb interconnect model for high-rate signal integrity”.  IEEE transactions on electromagnetic compatibility,  54(2), pp. 397–404.
  * [11] Tseng, F. H., Liang, T. T., Lee, C. H., Der Chou, L., and Chao, H. C., 2014.  “A star search algorithm for civil uav path planning with 3g communication”.  In 2014 Tenth International Conference on Intelligent Information Hiding and Multimedia Signal Processing, IEEE, pp. 942–945.
  * [12] Cui, W., Fan, J., Ren, Y., Shi, H., Drewniak, J. L., and DuBroff, R. E., 2003.  “Dc power-bus noise isolation with power-plane segmentation”.  IEEE transactions on electromagnetic compatibility,  45(2), pp. 436–443.
  * [13] Cheng, H.-D., Jiang, X. H., Sun, Y., and Wang, J., 2001.  “Color image segmentation: advances and prospects”.  Pattern recognition,  34(12), pp. 2259–2281.
  * [14] Pham, D. L., Xu, C., and Prince, J. L., 2000.  “Current methods in medical image segmentation”.  Annual review of biomedical engineering,  2(1), pp. 315–337.
  * [15] Haralick, R. M., and Shapiro, L. G., 1985.  “Image segmentation techniques”.  Computer vision, graphics, and image processing,  29(1), pp. 100–132.
  * [16] Wang, L., Fonseca, R., and Tian, Y., 2020.  “Learning search space partition for black-box optimization using monte carlo tree search”.  arXiv preprint arXiv:2007.00708.
  * [17] Buck, R. C., 1943.  “Partition of space”.  The American Mathematical Monthly,  50(9), pp. 541–544.
  * [18] Whitley, D., 1994.  “A genetic algorithm tutorial”.  Statistics and computing,  4(2), pp. 65–85.
  * [19] Jing, T., Lim, M. H., and Ong, Y. S., 2003.  “A parallel hybrid ga for combinatorial optimization using grid technology”.  In The 2003 Congress on Evolutionary Computation, 2003. CEC’03., Vol. 3, IEEE, pp. 1895–1902.
  * [20] Mühlenbein, H., 1992.  “Parallel genetic algorithms in combinatorial optimization”.  In Computer science and operations research. Elsevier, pp. 441–453.
  * [21] Lipowski, A., and Lipowska, D., 2012.  “Roulette-wheel selection via stochastic acceptance”.  Physica A: Statistical Mechanics and its Applications,  391(6), pp. 2193–2196.
  * [22] Hassanat, A., Alkafaween, E., Alnawaiseh, N., Abbadi, M., Alkasassbeh, M., and Alhasanat, M., 2016.  “Enhancing genetic algorithms using multi mutations: Experimental results on the travelling salesman problem”.  International Journal of Computer Science and Information Security,  14, 09, pp. 785–801.
  * [23] Gong, X.-Y., Su, H., Xu, D., Zhang, Z., Shen, F., and Yang, H.-B., 2018.  “An overview of contour detection approaches”.  International Journal of Automation and Computing,  15, 06, pp. 1–17.
  * [24] Huttenlocher, D. P., Klanderman, G. A., and Rucklidge, W. J., 1993.  “Comparing images using the hausdorff distance”.  IEEE Transactions on pattern analysis and machine intelligence,  15(9), pp. 850–863.
  * [25] Andoni, A., Indyk, P., and Krauthgamer, R., 2008.  “Earth mover distance over high-dimensional spaces.”.  In SODA, Vol. 8, pp. 343–352.
  * [26] Gelperin, D., 1977.  “On the optimality of a star”.  Artificial Intelligence,  8(1), pp. 69–76.
  * [27] Hu, J., and Sapatnekar, S. S., 2001.  “A survey on multi-net global routing for integrated circuits”.  Integration,  31(1), pp. 1–49.
  * [28] Ozdal, M. M., and Wong, M. D., 2007.  “Archer: a history-driven global routing algorithm”.  In 2007 IEEE/ACM International Conference on Computer-Aided Design, IEEE, pp. 488–495.
  * [29] Liao, H., Zhang, W., Dong, X., Poczos, B., Shimada, K., and Burak Kara, L., 2020\.  “A deep reinforcement learning approach for global routing”.  Journal of Mechanical Design,  142(6), p. 061701.
  * [30] Kruskal, J. B., 1956.  “On the shortest spanning subtree of a graph and the traveling salesman problem”.  Proceedings of the American Mathematical society,  7(1), pp. 48–50.
  * [31] Peterson, L. E., 2009.  “K-nearest neighbor”.  Scholarpedia,  4(2), p. 1883.
  * [32] Wayahdi, M. R., Ginting, S. H. N., and Syahputra, D., 2021.  “Greedy, a-star, and dijkstra’s algorithms in finding shortest path”.  International Journal of Advances in Data and Information Systems,  2(1), pp. 45–52.
  * [33] Dogar, M. R., Koval, M. C., Tallavajhula, A., and Srinivasa, S. S., 2014.  “Object search by manipulation”.  Autonomous Robots,  36(1), pp. 153–167.
  * [34] Dyer, C. R., 1980.  “Computing the euler number of an image from its quadtree”.  Computer graphics and image processing,  13(3), pp. 270–276.
  * [35] Abel, L. C., 1972.  “On the ordering of connections for automatic wire routing”.  IEEE Transactions on Computers,  100(11), pp. 1227–1233.
  * [36] Chen, H.-Y., and Chang, Y.-W., 2009.  “Global and detailed routing”.  In Electronic Design Automation. Elsevier, pp. 687–749.