electric hexadecapole moment $\Pi$, the single-particle nuclear shell model again predicts zero (similar to $Q$) because the valence nucleon is electrically neutral. We are not aware of any nuclear structure calculations for $\Pi$ of $^{173}\mathrm{Yb}$. We estimate $\Pi \approx Q^2 \approx 9 \, \mathrm{b^2}$ as both $Q$ and $\Pi$ arise due to nuclear deformation; we will take this value as fiducial in further computations.

From Eqs. (\[Eq:1st-order-ME\],\[eq:hpf-constants\]), we obtain the following first-order energy corrections, $$\begin{aligned}
    W_6^{(1)}&=&\frac{35}{4}A+\frac{1}{4}B+C+D  \, , \nonumber \\
    W_5^{(1)}&=&\frac{11}{4}A-\frac{37 }{140}B-\frac{109}{35}C-\frac{41}{7}D \, ,\nonumber  \\
    W_4^{(1)}&=&-\frac{9}{4}A-\frac{3}{10}B+\frac{46}{35}C+12 D\, , \nonumber \\
    W_3^{(1)}&=&-\frac{25}{4}A-\frac{1}{14}B+\frac{22}{7}C-\frac{44}{7}D \, , \label{Eq:W1-individual}\\
    W_2^{(1)}&=&-\frac{37}{4}A+\frac{1}{4}B+\frac{11}{35}C-11D \, ,\nonumber  \\
    W_1^{(1)}&=&-\frac{45}{4}A+\frac{15}{28}B-\frac{33}{7}C+\frac{99}{7}D \, \nonumber .\end{aligned}$$

The second-order corrections are computed from Eqs. (\[eq:2nd-order-general\],\[eq:2nd-order\]), where we keep magnetic dipole and electric quadrupole contributions. To streamline the notation, we introduce dipole-dipole, dipole-quadrupole, and quadrupole-quadrupole constants. These are defined for individual intermediate states $|\gamma'J'\rangle$,

$$\begin{aligned}
\label{eq:2nd-order-term}
    \eta_{\mu\mu}[\gamma'J']&=&\frac{(I+1)(2I+1)}{I}\frac{\mu^2{\ensuremath{\langle \gamma J|}} |T_1^e|{\ensuremath{|\gamma' {J{'}}\rangle}}^2}{E_{\gamma{J}}-E_{\gamma' J' }}\,,\nonumber \\
    \eta_{\mu{Q}}[\gamma'J']&=&{\frac{(I+1)(2I+1)}{I}}\sqrt{\frac{2I+3}{2I-1}}\times\frac{\mu{Q}{\ensuremath{\langle \gamma{J}|}}|T_1^e|{\ensuremath{|\gamma'{J'}\rangle}}{\ensuremath{\langle \gamma{J}|}}|T_2^e|{\ensuremath{|\gamma'{J'}\rangle}}}{E_{\gamma{J}}-E_{\gamma'{J'}}}\,, \\
    \eta_{QQ}[\gamma'J']&=&\frac{(2I+1)(I+1)(2I+3)}{4I(2I-1)}\frac{Q^2{\ensuremath{\langle \gamma{J}|}}|T_2^e|{\ensuremath{|\gamma' {J'}\rangle}}^2}
    {E_{\gamma{J}}-E_{\gamma' {J'}}}\,.\nonumber\end{aligned}$$

Eq. (\[eq:2nd-order\]) shows that we need to sum over all possible intermediate states obeying both the parity and the angular selection rules - that is, the parity of the ${\ensuremath{|\gamma{J}\rangle}}$ and ${\ensuremath{|\gamma{'}{J'}\rangle}}$ states has to be the same and $|J+J'|\geq{k}\geq|J-J'|$. Thus, for dipole-dipole 