Forbidden transitions in neutral and charged current interactions between low-energy neutrinos and Argon

TL;DR

This paper calculates low-energy neutrino interaction cross sections with argon, emphasizing the importance of forbidden transitions, using a HF+CRPA framework relevant for neutrino detection experiments.

Contribution

It introduces a detailed HF+CRPA model to accurately include forbidden transitions in neutrino-argon interactions at low energies.

Findings

Forbidden transitions significantly affect reaction strength.

Modeling must include forbidden transitions for accurate simulations.

Results are relevant for supernova neutrino detection experiments.

Abstract

Background: The study of low-energy neutrinos and their interactions with atomic nuclei is crucial to several open problems in physics, including the neutrino mass hierarchy, CP-violation, candidates of Beyond Standard Model physics and supernova dynamics. Examples of experiments include CAPTAIN at SNS as well as DUNE's planned detection program of supernova neutrinos. Purpose: We present cross section calculations for quasielastic charged current and neutral current neutrinos at low energies, with a focus on $^{40}$Ar. We also take a close look at pion decay-at-rest neutrino spectra, which are used in e.g. the SNS experiment at Oakridge. Method and results: We employ a Hartree Fock + Continuum Random Phase Approximations (HF+CRPA) framework, which allows us to model the responses and include the effects of long-range correlations. It is expected to provide a good framework to calculate forbidden transitions, whose contribution which we show to be non-negligible. Conclusions: Forbidden transitions can be expected to contribute sizeably to the reaction strength at typical low-energy kinematics, such as DAR neutrinos. Modeling and Monte Carlo simulations need to take all due care to account for the influence of their contributions.

Tables

Table 1. Table 1: The total charged–current ( ν e subscript 𝜈 𝑒 \nu_{e} , 56 Fe) cross section value, folded with a DAR electron neutrino spectrum, tabulated for various models.

	$⟨{\sigma }_{DAR}⟩$ (${10}^{-42}{\mathrm{cm}}^2$)
HF+CRPA	212.9
G–Matrix QRPA Cheoun et al. (2010)	173.5
Phenomenological Mintz (2002)	214
Hybrid Kolbe and Langanke (2001)	240
Hybrid Paar et al. (2011)	259
RHB+RQRPA Paar et al. (2011)	263
LFG+RPA Sajjad Athar et al. (2006)	277
QRPA Lazauskas and Volpe (2007)	352
Exp. (KARMEN) Kolbe et al. (1999a)	$256\pm 108\pm 43$

Table 2. Table 2: The total 56 Fe B ​ ( GT − ) B superscript GT \mathrm{B(GT^{-})} strength, tabulated for various models from Ref.  Paar et al. ( 2011 ) .

	$B\left({\mathrm{GT}}^{-}\right)$
HF+CRPA	8.8
GXPF1J	9.5
DD–ME2	11.3
SGII	12.3
SLy5	14.0
Exp.	9.9 $\pm$ 2.5

Equations

$$\frac{\mathrm{d}\sigma }{\mathrm{d}{T}_f\mathrm{d}{\Omega }_f}=$$

$$\times (v_{CC}{W}_{CC}+v_{CL}{W}_{CL}+v_{LL}{W}_{LL}$$

$$+v_T{W}_T\pm v_{T^{\prime }}{W}_{T^{\prime }}),$$

$${\mathcal{J}}_{\lambda }^{nucl}(\omega ,\mathbf{q})=⟨{\Phi }_{\text{f}}|{\hat{J}}_{\lambda }(\mathbf{q})|{\Phi }_{\text{0}}⟩.$$

$${\hat{J}}_{\lambda }(\mathbf{q})=\int \mathrm{d}\mathbf{x}{e}^{i\mathbf{x}\cdot \mathbf{q}}{\hat{J}}_{\lambda }(\mathbf{x}),$$

$$J_0(\mathbf{q})=$$

$$J_3(\mathbf{q})=$$

$$J_{\pm }(\mathbf{q})=$$

$${\hat{M}}_{J,0}^C=$$

$${\hat{T}}_{1,m}^E=$$

$$|{\Psi }_{RPA}^c⟩=\sum _{c^{\prime }}\left\{X_{c,c^{\prime }}|p^{\prime }{h}^{\prime -1}⟩-Y_{c,c^{\prime }}|h^{\prime }{p}^{\prime -1}⟩\right\},$$

$${\Pi }^{RPA}(x_1,x_2,E_{exc})={\Pi }^{(0)}(x_1,x_2,E_{exc})$$

$$+\frac{1}{\hslash }\int \mathrm{d}x\int \mathrm{d}{x}^{\prime }[{\Pi }^{(0)}(x_1,x,E_{exc})$$

$$\times \tilde{V}(x,x^{\prime }){\Pi }^{RPA}(x^{\prime },x_2,E_{exc})],$$

$$⟨p{h}^{-1}|\hat{O}|{\Phi }_{\text{0}}⟩\to v_h⟨p{h}^{-1}|\hat{O}|{\Phi }_{\text{0}}⟩,$$

$${\zeta }^2(Z^{\prime },E_f)=$$

$$\times$$

$$E_{eff}=E_f-V_c(0)=E\pm \frac{3}{2}\frac{Z^{\prime }\alpha \hslash c}{R},$$

$${\zeta }^2(Z^{\prime },E_f)=\frac{E_{eff}{k}_{eff}}{E_f{k}_f},$$

$$n(E_{{\nu }_e})=\frac{96E_{{\nu }_e}^2}{m_{\mu }^4}(m_{\mu }-2E_{{\nu }_e}),$$

$$n(E_{\mathrm{bar}{\nu }_{\mu }})=\frac{32E_{\mathrm{bar}{\nu }_{\mu }}^2}{m_{\mu }^4}\left(\frac{3}{2}{m}_{\mu }-2E_{\mathrm{bar}{\nu }_{\mu }}\right),$$

Full text

Purpose: We present cross section calculations for charged current and neutral current neutrino–nucleus scattering at low energies, with a focus on 40Ar. We also take a close look at pion decay–at–rest neutrino spectra, which are used in e.g. the SNS experiment at Oakridge.

Method and results: We employ a Hartree Fock + Continuum Random Phase Approximations (HF+CRPA) framework, which allows us to model the responses and include the effects of long–range correlations. It is expected to provide a good framework to calculate forbidden transitions, whose contribution we show to be non–negligible.

Conclusions: Forbidden transitions can be expected to contribute sizeably to the reaction strength at typical low–energy kinematics, such as DAR neutrinos. Modeling and Monte Carlo simulations need to take all due care to account for the influence of their contributions.

I Introduction

Within the realm of neutrino physics, one area of research that has a lot of potential for exciting new discoveries lies in the study of low-energy neutrinos. The interaction of these with atomic nuclei provides crucial information needed for several open questions in physics. DUNE DUN ; Acciarri et al. (2015) can, for example, as part of their low–energy (LE) neutrino program, distinguish between the two hierarchy scenarios in the detection of supernova neutrinos as explained in e.g. Ankowski et al. (2016), where mention is also made of the prospect to probe Beyond Standard Model (BSM) physics such as sterile neutrinos, neutrino magnetic moments and Non–Standard Interactions (NSI). Furthermore, the measurement of low–energy neutrinos is crucial to supernova research, where detecting the outgoing neutrinos can be used to study the underlying supernova dynamics and the interactions that the neutrinos engage in Ankowski et al. (2016). Finally, another avenue worth mentioning lies in the experimental measurement of coherent elastic neutrino-nucleus scattering (CE$\nu$NS) Scholberg (2015); coh ; Akimov et al. (2017).

Inelastic cross sections induced by low–energy neutrinos (10s of MeVs) have been measured in the past for carbon and iron by the LSND and KARMEN collaborations Auerbach et al. (2001); Maschuw (1998). Several experiments have been proposed for the near future to measure supernova neutrinos at higher precision, including the aforementioned DUNE as well as e.g. JUNO jun and Hyper-Kamiokande hyp . These require significant research and development efforts in order to properly calibrate the detectors that are to be used, including a better knowledge of the cross sections. One necessity is to have access to neutrinos of similar kinematics to those produced in supernovae. The Spallation Neutrino Source at Oakridge (SNS) Efremenko and Hix (2009) is an example of such a facility. As a byproduct of its primary purpose of producing neutrons, the SNS also creates neutrinos of several flavors. These are born out of pions decaying at rest (DAR), yielding monoenergetic muon neutrinos at $E_{\nu_{\mu}}=29.8$ MeV as well as electron neutrinos and muon antineutrinos, with their normalized flux pictured in Fig. 1, the well–known Michel spectra. Part of the CAPTAIN cap program is to have the CAPTAIN Liquid Argon Time Projection Chamber (LarTPC) detector run at the SNS Berns et al. (2013). This will yield crucial information on the neutrino-nucleus cross sections, but also on how to characterize these events in the analyses. A proper analysis of all the produced particles is a necessity to perform accurate calorimetry in SN neutrino experiments. Since it is of relevance to further discussions in this paper, we also make mention of the fact that specialized event generators are required for simulations, since most like NuWro nuw , NEUT Hayato (2009) and GENIE gen are tailored towards higher-energy scenarios. DUNE uses MARLEY mar , for example, to model charged current interactions between LE neutrinos and 40Ar.

Parallel to all of these developments lies a need to theoretically study and model cross-sections for charged current (CC) and neutral current (NC) interactions of SN neutrinos with nuclei. Past efforts to that effect include calculations performed in RPA frameworks such as Refs. McLaughlin (2004); Engel et al. (2003); Volpe et al. (2002); Auerbach et al. (1997); Suzuki and Sagawa (2003); Bandyopadhyay et al. (2017); Cheoun et al. (2011), shell model Hayes and Towner (2000); Suzuki et al. (2006); Suzuki and Kajino (2013); Kostensalo et al. (2018) (or hybrid models of the two Kolbe and Langanke (2001); Kolbe et al. (1999a); Suzuki et al. (2009); Suzuki and Honma (2013)) as well as QRPA Paar et al. (2008); Samana and Bertulani (2008); Cheoun et al. (2010), relativistic (R)QRPA Paar et al. (2011) and local Fermi Gas–based RPA approaches Nieves et al. (2004); Sajjad Athar et al. (2006); Singh et al. (1998). We also mention previous CRPA results in Kolbe et al. (1999b); Jachowicz et al. (2002a); Jachowicz and Heyde (2003). Comparisons between shell model, RPA and CRPA results were performed in Volpe et al. (2000); Kolbe et al. (2003). Phenomenological calculations of cross sections were also provided, in Refs. Pourkaviani and Mintz (1990); Fukugita et al. (1988); Mintz (2002). In this paper, we present cross section calculations for low–energy CC and NC neutrino–nucleus scattering in a HF+CRPA approach, with a focus on 40Ar.

II Model

We will briefly discuss the model used for our calculations. The cross section for inclusive electroweak scattering of neutrinos off atomic nuclei is given by the familiar expression:

[TABLE]

where $\sigma_{X}$ is $\left(\frac{G_{F}\cos{\theta_{c}}}{2\pi}\right)^{2}$ for CC scattering and $\left(\frac{G_{F}}{2\pi}\right)^{2}$ for NC. In turn, $G_{F}$ is the Fermi constant that encodes the strength of the weak interaction, while $\cos{\theta_{c}}$ is the Cabibbo angle. The factor $\zeta^{2}(Z^{\prime},E_{f})$ is related to the Coulomb interaction between the escaping lepton and the residual nucleus. In the case of NC, it is simply 1. We will come back to the specifics of this factor further on. The cross section also differs for neutrinos and antineutrinos due to the V–A nature of the electroweak interaction, as shown by the $\pm$–sign in the transverse interference term. The $v$–factors are the leptonic functions, while the $W$–factors are the nuclear response functions, which contain all the nuclear information on the reaction to the electroweak probe. They are a function of the transition amplitude between the initial $|\Phi_{\textrm{0}}\rangle$ and final $|\Phi_{\textrm{f}}\rangle$ state:

[TABLE]

This nuclear current is the Fourier transform of the current in coordinate space:

[TABLE]

where $\omega$ and $\boldsymbol{q}$ are the energy and momentum transfer, respectively. One performs a multipole expansion of the nuclear current, which can be used to exploit angular momentum algebra. If one chooses the z–direction along the momentum transfer $\boldsymbol{q}$:

[TABLE]

Expressions for the nuclear one–body current and derived Coulomb ($\hat{M}_{J,L}^{C}$), longitudinal ($\hat{L}_{J,L}^{L}$), electric ($\hat{T}_{J,L}^{E}$) and magnetic ($\hat{T}_{J,L}^{M}$) multipole operators can be found in e.g. Walecka (2004). Worth mentioning at this stage in the context of low–energy neutrinos, is the ’allowed approximation’. Under this approximation, valid at low energies and typically used in beta–decay calculations, it is customary to employ both the long–wavelength limit $q\rightarrow 0$, and to assume that one is dealing with slow nucleons $p_{N}/m_{N}\rightarrow 0$. Under these assumptions, one can show that the only terms surviving in the multipole decomposition of the nuclear current for CC interactions (similar considerations hold for the NC current) are the following:

[TABLE]

where $F_{1}$ and $G_{A}$ are the Fermi and axial form factors, respectively. These are the well–known Fermi and Gamow–Teller transition operators. At sufficiently low energies, they provide an adequate description of the nucleus’ response to an electroweak probe. Because they are responsible for the largest part of the reaction strength, the transitions they give rise to are known as ’allowed’, with the $J^{P}$ quantum numbers of the operators $0^{+}$ and $1^{+}$, respectively. Conversely, higher–order transitions are often referred to as ’forbidden’ transitions. As we will show, the latter are still important at the energies considered in this context. In this paper, the nuclear responses are calculated in a Continuum Random Phase Approximation (CRPA), where long–range correlations and collective excitations of the nucleus are taken into account. We mention that in these calculations, we use the free-nucleon value for the axial coupling of $g_{A}=1.27$, whereas some models use an effective quenched value of around $g_{A}=1.00$. This issue is discussed in Ref. Pastore et al. (2018), where it is argued that the need for such a quenching may appear as a consequence of an absence of nuclear correlations or an insufficiently large model space in the modeling.

Within the HF+CRPA approach, excited states of the nucleus are described as coherent superpositions of particle–hole and hole–particle excitations out of the correlated ground state:

[TABLE]

where c contains all the quantum numbers to unambiguously label an excitation channel. The local RPA polarization propagator, containing all the information on excited states, is implicitly defined as

[TABLE]

with $\tilde{V}(x,x^{\prime})$ the antisymmetrized residual interaction, $E_{exc}$ the excitation energy of the nucleus and $x$ a shorthand for spatial, spin and isospin coordinates. The residual interaction used is the same Skyrme force that is employed to calculate the single–particle wave functions, so that the scheme is self–consistent. Partially filled shells such as those present in 40Ar are dealt with by including occupation probabilities in the transition amplitudes:

[TABLE]

with $\hat{O}$ representing a general one–body operator and $v_{h}^{2}$ being the occupation probability of the shell $h$. Further details and previous applications of this model can be found in Refs. Ryckebusch et al. (1988, 1989); Jachowicz et al. (1999); Jachowicz et al. (2002b, a); Jachowicz et al. (2004); Jachowicz and McLaughlin (2006); Jachowicz et al. (2008); Pandey et al. (2015a, 2016); Van Dessel et al. (2018); Nikolakopoulos et al. (2018). While our calculations predominantly pertain to describing reactions in the Giant Resonance (GR) region due to the low neutrino energies, the HF+CRPA model described above is also successful at describing the quasielastic region at higher energies, such as those relevant to baseline experiments Pandey et al. (2014, 2015a, 2015b, 2016); Van Dessel et al. (2018).

Finally, we turn our attention to the topic of Coulomb interactions between the outgoing lepton and the residual nucleus in charged–current interactions. At low energies, this is expected to have a strong effect on the cross section due to the low momentum of the charged outgoing lepton. In our discussion we will consider two effective schemes to account for the Coulomb distortion. At low outgoing lepton energies, of the order of magnitude applicable to e.g. beta decay, one can make use of the Fermi function Engel (1998):

[TABLE]

where $R\approx 1.2A^{1/3}\mathrm{fm}$ is the nuclear radius, $\gamma_{0}=\sqrt{1-(\alpha Z^{\prime})^{2}}$, $E_{f}$ is the outgoing lepton’s energy, $p_{f}$ outgoing momentum and $\eta=\pm\frac{\alpha Z^{\prime}c}{v}$. with $+$ and $-$ for neutrinos and antineutrinos respectively. Similarly, the final nuclear charge $Z^{\prime}$ is equal to $Z+1$ or $Z-1$ for $\nu$/$\mathrm{bar}{\nu}$, respectively. This approximation assumes that the outgoing leptons only contribute sizeably as an s–wave to the reaction strength, and is therefore not applicable at higher outgoing lepton energies Engel (1998). Therefore, in these regimes, we will consider a different scheme, the modified effective momentum approximation (MEMA), detailed in Ref. Engel (1998), where the energy and momentum of the final lepton is shifted to an effective value by the Coulomb energy in the center of the nucleus:

[TABLE]

which introduces a factor in the cross section that accounts for the change in phase space:

[TABLE]

and induces a shift in the momentum transfer $q\rightarrow q_{eff}$ in the calculation of the amplitudes in Eq. (2). In practice, we will interpolate between the two schemes. This approach consists of taking for each value of $\omega$ in the differential cross section the smallest or highest value of the two, respectively for neutrinos and antineutrinos where the effect of the Coulomb interaction is to increase resp. decrease both the differential and the integrated cross section Engel (1998). In short: we take the value of $\zeta^{2}(Z^{\prime},E_{f})$ that is closest to unity. The effect of this is shown in Fig. 3. At low incoming neutrino energies, the interpolation scheme matches that of the Fermi function, whilst at higher energies, where reaction strength contains increasingly more contributions from events with high momentum of the final lepton, the scheme matches that of the MEMA. We show this for several nuclei. For 12C, it is clear that the Fermi function is a good approximation over the whole energy range, while for 40Ar and 56Fe, the Fermi function and MEMA schemes cross between 60 and 70 MeV. In general, the heavier the nucleus, the more limited the Fermi function will be in its area of applicability. A similar comparison of Coulomb schemes was performed in Ref. Sajjad Athar et al. (2006).

III Results

Before we cover the results for neutrino–argon calculations, we wish to compare our model predictions with other recent studies at these energies for various nuclei. In Fig. 4 we confront the HF+CRPA results for 56Fe to those of other models. The comparison is fine, except that the strength lies comparatively low for lower energies. This is because of the fact that, within the CRPA, discrete excitations are not taken into account. We can also take a look at the cross section predictions for DAR kinematics, by folding this cross section with the appropriate neutrino spectrum. We compare this to the KARMEN experimental results and various other theoretical models in Table 1. Our result is compatible with experimental results, but has a lower predicted value than most other theoretical models. This is due to the fact that, once again, discrete excitations are not included in the reaction strength. Finally, we compare the $\mathrm{B(GT^{-})}$ strength with other models discussed in Paar et al. (2011). Firstly, this strength, as a function of the nuclear excitation $E_{x}$ is shown in Fig.  5. The HF+CRPA results lack the discrete excitation strength, but reproduces the large peak of the QRPA results relatively well. Secondly, we similarly compare the total strength in Table 2. HF+CRPA underpredicts the other models, likely once again due to the absence of discrete excitations, but still falls within the broad experimental error interval.

In Fig. 6 we compare with other model predictions for argon. A similar comparison holds. The results from MARLEY match surprisingly well with those from Sajjad Athar et al. (2006), in spite of using a different, Fermi Gas–based RPA model. The other two models predict lower strength at higher energies. For Ref. Kostensalo et al. (2018) in particular, using a shell model, it is argued that this is due to a limited valence space for excited states at higher energies.

We now provide HF+CRPA model predictions for neutrino–40Ar cross sections at low incoming neutrino energies. We first present total cross sections as a function of incoming energy. In figures 7, 8, 9 and 10 we show the reaction strength for CC and NC, for both incoming $\nu_{e}$ and $\mathrm{bar}{\nu}_{e}$.

In the plots the total cross sections are shown for different multipole cutoffs $J$. We also show the contribution of natural ($0^{+},1^{-},2^{+},3^{-}$, …) and unnatural ($0^{-},1^{+},2^{-},3^{+}$, …) parity transitions separately. In comparing the case of CC and NC for neutrinos and antineutrinos, other than the overall magnitudes of the cross sections, the qualitative features concerning the relative contributions from different multipolarities are similar. The $0^{+}$ and $1^{+}$ transitions correspond with the allowed approximation. One can appreciate the non–trivial contributions arising from the forbidden $1^{-}$ and $2^{-}$ transitions in the continuum channels. The importance of these two in particular was similarly noted in Ref. Kostensalo et al. (2018). Even higher order multipoles start becoming visible above 50 MeV. The $0^{+}$ reaction strength is negligible for excitations into the continuum. In comparing the CC and NC result, one notices that the difference in strength between neutrinos and antineutrinos is noticeably larger for the former. A key reason for this is the Coulomb interaction with the outgoing lepton, which enhances or decreases the cross section for neutrinos and antineutrinos, respectively.

Next, we cover the energy dependence of the contributions of different multipole moments. We first take a look at the differential cross section $\frac{\mathrm{d}\sigma}{\mathrm{d}\omega}$ for 40Ar in Fig. 11 for several incoming energies. For clarity, we do not show the $0^{-}$ and $0^{+}$ contributions, since these are small. At 25 MeV incoming energy, we can see that the cross section is quite well described in an allowed approximation, with the differential cross section almost coinciding with the $1^{+}$ transition strength, except for a small difference in the tail. At 50 MeV, this is no longer an adequate approximation, with significant strength coming from the $1^{-}$ and $2^{-}$ transitions, affecting the shape of the excitation spectrum as well. Finally, at 75 MeV, even higher order multipolarities become non-negligible. In Figs. 12 and 13 we investigate the A–dependence of the cross sections. Comparing Figs. 7 to 12 and 13, one observes that even for a relatively small nucleus as 12C, forbidden transitions are non–negligible at higher energies, but not to nearly the same degree as they are for 40Ar and 56Fe, both of which show a qualitatively similar picture, not differing too significantly in size. In general, the larger the target nucleus, the more multipoles will be necessary for convergence of the cross section.

In order to appreciate the effect that the individual multipole contributions will have in an experimental context, it is instructive to take a look at folded cross sections. For now, we make use of the spectra displayed in Fig. 1, yielded by $\pi^{+}$ pions decaying at rest (DAR), with associated fluxes

[TABLE]

and

[TABLE]

with $E_{\nu}\in\left[0,m_{\mu}/2\right]$ and $m_{\mu}\approx$ 105.7 MeV the rest mass of a muon. We fold the cross sections on display in figures 7, 8, 9 and 10 (except for CC ($\mathrm{bar}{\nu}_{\mu},{}^{40}$Ar), which is not kinematically available) with the appropriate fluxes of Eqs. 12 and 13. The results, plotted individually for each multipole, are shown in Figs. 14, 15 and 16. The plots demonstrate that for argon nuclei, forbidden transitions can be expected to contribute sizeably to the reaction strength in the energy regimes of DAR neutrinos. All three paint a qualitatively similar picture in terms of relative contributions, except for the case of NC ($\mathrm{bar}{\nu}_{\mu},{}^{40}$Ar), where the $\mathrm{bar}{\nu}_{\mu}$ flux peaks at a higher value, as seen in Fig. 1. Therefore the $2^{-}$ transition strength has a larger relative weight. And although small, even the $2^{+}$, $3^{-}$ and $3^{+}$ transitions cannot be neglected when high precision is desired. We conclude that working strictly in an allowed approximation is not desirable, and forbidden transitions must be properly accounted for in the modeling of low–energetic neutrino–nucleus scattering events. Therefore it is worth noting that the Monte Carlo generator MARLEY mar employs Fermi and Gamow–Teller transitions in its modeling of charged–current neutrino–40Ar scattering, but does not go beyond this Gardiner (2017). In light of our results, the implementation of forbidden transitions should be considered important.

IV Summary

Argon represents an important nuclear scattering target in low–energy neutrino research. We have presented calculations for ($\nu,{}^{40}$Ar) scattering. We employed a Hartree Fock + Continuum Random Phase Approximations (HF+CRPA) framework, which allows us to model the responses and include the effects of long–range correlations. The Coulomb interaction of the outgoing lepton is accounted for in an effective way by interpolating between the Fermi function, valid at low outgoing charged lepton momenta, and the Modified Effective Momentum Approach (MEMA) at higher ones. Subsequently comparing model results for 56Fe nuclei, one notices that the model predictions compare favorably. Finally, we presented our predictions for 40Ar, for both neutral and charged current events, where the important role played by forbidden transitions is evident. Even in a typical experimental scenario of a $\pi^{+}$ decay–at–rest neutrino spectrum, this is still the case, with strong contributions from the $1^{-}$ and $2^{-}$ multipoles, and smaller contributions for higher multipolarities.

Acknowledgements.

This work was supported by the Research Foundation Flanders (FWO–Flanders). The computational resources (Stevin Supercomputer Infrastructure) and services used in this work were provided by the VSC (Flemish Supercomputer Center), funded by Ghent University, FWO and the Flemish Government – department EWI.