Control defeseance by anti-alignment in the excited state

TL;DR

This paper predicts that in the excited state of certain molecules, strong fields induce anti-alignment due to opposite polarizability signs, challenging existing laser-control strategies.

Contribution

It reveals a fundamental difference in polarizability effects between ground and excited states, showing anti-alignment dynamics in excited molecules under strong fields.

Findings

Anti-alignment occurs in excited states of H2+ and similar molecules.

Polarizability signs invert from ground to excited states, causing opposite alignment behavior.

Laser control strategies effective in ground states may fail in excited states.

Abstract

We predict anti-alignment dynamics in the excited state of H2+ or related homonuclear dimers in the presence of a strong field. This effect is a general indirect outcome of the strong transition dipole and large polarizabilities typically used to control or to induce alignment in the ground state. In the excited state, however, the polarizabilities have the opposite sign than in the ground state, generating a torque that aligns the molecule perpendicular to the field, deeming any laser-control strategy impossible.

Equations

$$i\frac{\partial }{\partial t}\psi (\tilde{\mathbf{x}})=-\frac{1}{2}\sum _j\frac{{\partial }^2}{\partial {\tilde{x}}_j^2}\psi (\tilde{\mathbf{x}})+\left(V_{sc}+V_{int}\right)\psi (\tilde{\mathbf{x}})$$

$$\frac{d^2}{d{t}^2}{\tilde{X}}_j^{(\alpha )}=-\frac{1}{M_{\alpha }}⟨\psi (\tilde{\mathbf{x}})\left|\frac{\partial (V_{sc}+V_{int})}{\partial {\tilde{X}}_j^{(\alpha )}}\right|\psi (\tilde{\mathbf{x}})⟩$$

$$V_{sc}=-\sum _{\alpha }{Z}_{\alpha }{\left[\sum _j{\left({\tilde{x}}_j-{\tilde{X}}_j^{(\alpha )}\right)}^2+{ϵ}^2\right]}^{-\frac{1}{2}}+Z_1{Z}_2{\left[\sum _j{\left({\tilde{X}}_j^{(2)}-{\tilde{X}}_j^{(1)}\right)}^2\right]}^{-1}$$

$$V_{int}=\sum _j{E}_j(t)\left({\tilde{x}}_j-Z_2{\tilde{X}}_j^{(2)}-Z_1{\tilde{X}}_j^{(1)}\right)$$

$$V_{sc}=-{\left[(x+\frac{r}{2}\cos \theta {)}^2+(y+\frac{r}{2}\sin \theta {)}^2+{ϵ}^2\right]}^{-1/2}-{\left[(x-\frac{r}{2}\cos \theta {)}^2+(y-\frac{r}{2}\sin \theta {)}^2+{ϵ}^2\right]}^{-1/2}+\frac{1}{r}$$

$$\frac{d^2}{d{t}^2}r=\frac{1}{M}\left(V_1^3-V_2^3\right)\left[⟨x⟩\cos \theta +⟨y⟩\sin \theta \right]+\frac{r}{2M}\left(V_1^3+V_2^3\right)-\frac{2}{M{r}^2}$$

$$\frac{d^2}{d{t}^2}\theta =\frac{1}{Mr}\left(V_1^3-V_2^3\right)\left[⟨y⟩\cos \theta -⟨x⟩\sin \theta \right]$$

$$\frac{d^2}{d{t}^2}{\mathbf{X}}^{CM}=\frac{1}{4M}\mathbf{E}(t)$$

$$i\frac{\partial }{\partial t}\psi (\mathbf{x})=-\frac{1}{2}\sum _j\frac{{\partial }^2}{\partial x_j^2}\psi (\mathbf{x})+\left[V_{sc}(\mathbf{x},r(t),\theta (t))+\mathbf{E}(t)\mathbf{x}-\mathbf{E}(t){\mathbf{X}}^{CM}(t)\right]\psi (\mathbf{x})$$

$$V_{\mathrm{ind}}=-\frac{1}{4}{E}^2(t)\left[{\cos }^2\theta \left({\alpha }_{\parallel }-{\alpha }_{\perp }\right)+{\alpha }_{\perp }\right],$$

$${\alpha }_{\parallel }^{(1)}=-2\sum _{n>1}\frac{⟨{\psi }_1|{\mu }_x|{\psi }_n⟩⟨{\psi }_n|{\mu }_x|{\psi }_1⟩}{V_1(r)-V_2(r)}\approx \frac{r^2}{2\Delta (r)}>0$$

Full text

I Introduction

Several experiments have reported successful control on the yields and products of photodissociation or other photo-induced processes using short and strong laser pulses Assion et al. (1998); Levis et al. (2001); Daniel et al. (2003); Brixner et al. (2004); Nuernberger et al. (2010). In principle, strong non-resonant fields can cleanly (i.e. adiabatically) manipulate the transition from reactants to products by distorting the reaction landscape through dynamical Stark-shifts that are sensitive to the reaction coordinateSussman et al. (2006); Kim et al. (2012); Corrales et al. (2014).

These manipulations are most effective when the electronic forces and gradients are small or comparable to the laser-induced couplings, allowing to shape true light-induced potentials, commonly referred to as LIPs Yuan and George (1978); Bandrauk and Sink (1981). But the limiting factor is the presence of the coupling, which in the dipole approximation is of the form $\mu_{ij}(r)E(t)\cos\theta$, where $E(t)$ is the laser amplitude, $\mu_{ij}(r)$ the transition dipole, and $\theta$ the angle between the molecular axis and the field polarization.

From the electronic structure point of view, effective quantum control relies on finding accessible pairs of states where the transition dipole does not decay to zero along the reaction coordinateSola et al. (2018). This is typically the bottleneck of many control scenarios. With limited pulse bandwidth (few lasers operate at enough power that most photon wavelengths are within some narrow ranges of frequencies) the wave packet is mostly affected by one or few transition dipoles to guide its motion, whereas the potential energy surface is quagmired by conical intersections where the dipole is essentially zeroYarkony (1996); Worth and Cederbaum (2004); Domcke et al. (2011).

In simple molecules (i.e. diatomic molecules) one can mostly avoid these problems. Indeed, aligned homonuclear diatomic cations provide the ultimate examples of highly-controllable moleculesChang et al. (2016). The symmetry-induced charge-transfer at dissociation is responsible for a huge transition dipole between the pair of states that instead of decaying asymptotically as the states become degenerate in the reaction channel, it increases linearly with the internuclear distance. This dipole can then be used to bind the molecule at any bond length, in a process known as LAMB (laser-adiabatic manipulation of the bond)Chang et al. (2003, 2004); Sola (2004); Gonzalez-Vazquez et al. (2006); Chang et al. (2010, 2011).

There is an additional source acting against laser control, not due to the vibration but to the rotation and orientation of the molecule with respect to the field, as noted in the $\cos\theta$ factor of the interaction. Within the LIPs framework, at $\theta=\pi/2$ the LIPs show light-induced conical intersections (LICIs) that often provoke undesired transitionsMoiseyev et al. (2008); Halász et al. (2012); Demekhin and Cederbaum (2013); Halász et al. (2015). In addition, near this point only molecular electronic gradients act on the nuclear wave packet, hence dissociative states lose their laser-induced binding properties. However, it is well known that lasers (and specially strong non-resonant lasers) force the alignment of the molecular axis with the polarization axis. Even when no prior alignment is pursued or measured, it is usually believed that the concomitant alignment is helping in the control of the dynamics or, at least, is a resource that the controller may exploit Balint-Kurti et al. (2008); Brif et al. (2012).

But while this is the case in the ground state, where one can easily show by second order perturbation theory that the static polarizabilities in the ground state have to be negative, inducing alignment, it is not necessarily so in the excited states. In this work we will show that precisely in diatomic homonuclear cations the polarizabilities in the excited state are basically just the opposite of the polarizabilities in the ground state, inducing anti-alignment. That is, the effective potential for the rotation in the presence of the field has a maximum when the molecular axis is aligned with the transition dipole and a minimum when it is perpendicular to the field. Then basically all the molecules that are not initially aligned with the field will dissociate in the excited state, defeating any control strategy.

The paper relies on a $2$-dimensional model Hamiltonian with soft-core Coulomb potentialsJavanainen et al. (1988); Su and Eberly (1991); Kulander et al. (1996), where the nuclear motion is included using the Ehrenfest approachTully (1998, 2012). For the aligned molecules (in a effectively $1$-dimensional treatment) the model gives results that are in agreement with those of full quantum (electron and nuclear) calculations with the same HamiltonianChang et al. (2019). But we expect, as explained in the paper, the validity of the results to be far more general, as well as the process of anti-alignment.

In the remaining of this work we first introduce the model Hamiltonian. Then we analyze the results obtained in the ground state of molecules at random orientations and in the excited state of molecules initially aligned with the field, before studying the dynamics of molecules initially misaligned with the field.

II Model

Our goal is to design a simple model that can provide qualitative predictions of quantum control in one-electron systems, treating strong field laser couplings with polarized lasers in non-aligned diatomic (or planar) molecules, non-adiabatic couplings and ionization on equal footing. To that end we use soft-core Coulomb potentials in a $2$-D grid, where the nuclear motion is incorporated in the Ehrenfest approach. Calling $\tilde{X}_{j}^{(\alpha)}$ ($j=1,2$ for the $X$, $Y$ coordinates, $\alpha=1,2$ for the nuclei of mass $M_{\alpha}$ and charge $Z_{\alpha}$) the Cartesian coordinates of the nuclei in the plane of motion, and equivalently $\tilde{x}_{j}$ the electron coordinates, collectively written as $\tilde{\bf{X}}$ and $\tilde{\bf{x}}$ respectively, and $\psi(\tilde{{\bf x}})$ the electron wave function, the electronic degrees of freedom obey the time-dependent Schrödinger equation (in atomic units)

[TABLE]

while for the nuclear degrees of freedom we use the Hellmann-Feynman force approximationTully (1998, 2012),

[TABLE]

where the soft-core Coulomb potential,

[TABLE]

has analytic derivatives. We have chosen $\epsilon=1/\sqrt{2}$. In the presence of a linearly polarized external field, ${\bf E}={\bf i}E_{1}+{\bf j}E_{2}$, in the dipole approximation, the interaction potential is given by

[TABLE]

In this work both nuclei will be protons ($Z_{\alpha}=1$, $M_{\alpha}=M$). It is convenient to decouple the internal and center of mass motion and work with coordinates relative to the center of mass, not only to reduce the number of variables, but mainly to avoid having to work with very large grids for the electronic coordinates, as the whole molecular cation moves in the gradient of the field. In order to do so we first disregard the effect of the magnetic field. Calling the nuclear center of mass coordinates ${\bf X}^{CM}$ and defining the distances with respect to ${\bf X}^{CM}$, ${\bf x}=\tilde{\bf x}-{\bf X}^{CM}$, ${\bf X}=\tilde{\bf X}-{\bf X}^{CM}$, we obtain the same equation for $V_{sc}({\bf x},{\bf X})$ as before [Eq.(4)], removing the tilde of the variables. If we neglect the difference between the nuclear center of mass and the molecule’s center of mass, then $V_{int}={\bf E}(t)\left({\bf x}-{\bf X}^{CM}\right)$. Within the same level of approximation, we also use the electron mass for the reduced mass of the system and neglect the mass-polarization term. Hence, the Ehrenfest equation for $X_{j}^{(\alpha)}$ remains the same as previously [Eq.(3)], but dropping the tilde. Now the equations depend only on the gradient of $V_{sc}$, while those for ${\bf X}^{CM}$ depend only on the gradient of $V_{int}$. For instance, using polar variables, $X_{1}^{(1)}=-r/2\cos\theta$, $X_{2}^{(1)}=-r/2\sin\theta$ ($X_{j}^{(2)}=-X_{j}^{(1)}$),

[TABLE]

$\equiv V_{1}+V_{2}+r^{-1}$. So,

[TABLE]

[TABLE]

where $\langle x\rangle\equiv\langle\psi({\bf x})|x|\psi({\bf x})\rangle$ and $\langle y\rangle\equiv\langle\psi({\bf x})|y|\psi({\bf x})\rangle$ are the average electronic positions, obtained by integrating the electronic wave function in the grid. While for the nuclear center of mass coordinates,

[TABLE]

which can be immediately integrated. In the new variables, the TDSE becomes

[TABLE]

The last term in the parenthesis is a time-dependent scalar that can be taken away from the TDSE by a unitary transformation of the wave function, which receives a time-dependent phase due to the accelerated motion of the center of mass of the charged molecule in the presence of the field, $\varphi(t)=\int_{0}^{t}{\bf E}(t^{\prime}){\bf X}^{CM}(t^{\prime})dt^{\prime}$, where ${\bf X}^{CM}(t)$ is obtained from the solution of Eq.(9).

In our simulations, the initial electron wave function is obtained by imaginary time propagation. To obtain the ground state for a given nuclear arrangement we initially propagate a Gaussian wave function. To obtain the $n$ excited state we start with a Gaussian times $x^{n}$ and filter all the lower-energy eigenstates previously obtained. For convenience, the molecule is initially aligned in the $x$ axis of the grid ($X_{2}^{(1)}=X_{2}^{(2)}=0$) while the laser polarization is appropriately rotated by an angle equal to $-\theta$. The dynamics is solved in a square grid of $256\times 256$ points between $-35$a0 and $35$a0 using the second-order split-operator method. The kinetic energy operator is evaluated by fast-Fourier transform and we use imaginary absorbing potentials at the grid edgesMacías et al. (1994).

To represent the behavior of an ensemble of molecules, we typically perform $200$ simulations starting with different initial conditions obtained from random sampling of the Wigner distribution of a Gaussian wave packet. This is either the ground state of the molecule or a wave packet in the excited state chosen at a particular position and with a particular width that resembles the wave packet obtained after ionizing H2 and exciting it with the third harmonic of an $800$ nm laser Chang et al. (2013a, 2015). In all the results of this work, we assume that all molecules are initially non-rotating, that is, the total angular momentum is zero.

III Results

In this work we will analyze the orientation, vibration and electron dynamics (the generation of dipoles) in H${}_{2}^{+}$ in the presence of a strong static field. Although the numerical model assumes some approximations, we expect the physical picture regarding the bound states dynamics to emerge in qualitative agreement with that of a fully quantum treatment. In particular, the average results obtained from an ensemble of trajectories (weighted by a Wigner distribution) using the Ehrenfest approximation for a $1$-D model of H${}_{2}^{+}$ with a soft-core Coulomb potential are very similar to those obtained by solving the quantum $1+1$D (electron-nuclear) TDSE, both when we start in a superposition of electronic states in the absence of the field, and when we start in the excited state (or the dressed electronic state) in the presence of the fieldChang et al. (2019).

The $2$-D electronic model has new features: it enables to calculate (and observe the influence of) states with orbital angular momentum, although without the right degeneracy, and more importantly, it allows to study the role of orientation. The predictions of the model for the dynamics in the ground state, in the presence of a strong static field (or a low-frequency continuous-wave laser) are consistent with the well-known process of alignment. In Figure 1(a) we show the dynamics of the average bond length $\langle r(t)\rangle$ and in Figure 1(b) the average orientation of the molecular axis with respect to the field, $\langle\theta(t)\rangle$, for a field amplitude of $E_{0}=0.02$ a.u. In a particular trajectory (orange line), the molecular axis is initially $\pi/4$ degrees with respect to the electric-field polarization (in the $x$ axis), while the initial bond length is $2.4$ a0. We also show the average results for an ensemble of $200$ molecules with random orientations with respect to the field and a Wigner distribution of the internal nuclear coordinates that corresponds to a Gaussian wave packet centered at $r_{0}=2.4$ a0 with width $\sigma=0.5$ a0, close to the ground state equilibrium bond distance for the $2$-D soft-core coulomb potential model. While in a particular trajectory the vibration is almost decoupled from the rotation in the ground state, given that the internuclear distance does not vary greatly and the moment of inertia of the molecule remains fairly constant, the averaged internuclear distances quickly (after the first oscillations) attenuate the vibrational motion, which is induced by the anharmonicity of the potential, and hence is sensitive to molecular orientation. There is also a small additional effect caused by bond-softening in the ground state due to the strong fieldBucksbaum et al. (1990); Allendorf and Szöke (1991); Jolicard and Atabek (1992); Giusti-Suzor et al. (1995). A single trajectory exhibits field-induced alignment, which shows in the pendular dynamics. The ensemble of molecules remains uniformly distributed over the orientations at all times.

When the molecule is in the first excited electronic state in the presence of the field, there is bond-hardeningGiusti-Suzor and Mies (1992); Yao and Chu (1992); Zavriyev et al. (1993); Aubanel et al. (1993); Chang et al. (2013b). The electron motion is correlated to the vibrational motion such that the electron moves away with the proton against the field, reaching a classical turning point and returning to the initial position. This was shown in a fully-quantum $(1+1)$-D model for a molecule aligned with the fieldChang et al. (2013a, 2015) and reproduced with a 1D Ehrenfest model of the electron moving in the molecular axisChang et al. (2019). Figure 2 depicts the results using the $2D$-Ehrenfest model for an ensemble of $200$ molecules initially aligned with a field of $E_{0}=0.02$ a.u. The initial nuclear coordinates were obtained from the Wigner distribution of a Gaussian state centered at $r=4$ a0, with different widths, corresponding to the conditions tested in the full quantum calculation of Chang et al. (2015). For perfect initial alignment, $\theta(t)=0$ and the molecules remain aligned for all times. As a particular trajectory, we show $r(t)$ starting at $r(0)=4$ a0 assuming zero initial vibrational kinetic energy for the nuclei. The average $\langle r(t)\rangle$ reproduces a decay in the oscillations, due to dephasing. The amplitude of the bond dynamics depends on the choice of the initial Wigner distribution. In our example, the narrower initial distribution ($\sigma=a_{0}/2$) reaches larger bond distances. This is because it corresponds to a broader initial momentum distributionChang et al. (2019). For the same reason the dephasing is somehow stronger for this distribution and the decay of the oscillation is larger than for the $\sigma=a_{0}/\sqrt{2}$ case. In spite of the dephasing, the first oscillations are clearly visible marking the regime where the vibrational dynamics remains coherent for the whole ensemble of molecules, in agreement with the full quantum calculation.

However, the dynamics differs completely when the molecule is not initially aligned with the field in the excited state. In this case, the main driver of the dynamics is the torque against the field that is exerted on the molecular axis. This anti-alignment weakens the coupling and eventually leads to dissociation when the molecular axis is perpendicular to the field. In figure 3 we show several examples of the dynamics for different initial orientations ($\theta(0)=\pi/4,\pi/10$ and $\pi/40$) when the pulse amplitude is $E_{0}=0.02$ a.u., $\langle r(0)\rangle=4$a0 and the nuclei start with zero vibrational kinetic energy in the excited electronic state.

The vibrational period in the excited LIP is approximately $50$ fs (see Fig.2). While the bond length enlarges initially as expected (see Fig.3a), the molecule is rotating against the field (see Fig.3b), making $\theta$ larger until there is no coupling at all, at $\theta=\pi/2$, and the molecule dissociates. As the coupling with the field becomes weaker, the LIP is becoming flatter and the bond recovery is nonexistent when the molecular axis is initially very misaligned with the field (e.g. $\theta(0)=\pi/4$), causing a fairly direct photodissociation. For small initial misalignment ($\theta(0)=\pi/10$ or $\pi/40$) the coupling is stronger and one can observe some (partial) vibration in the LIP. In addition, if the bond is not too stretched when $\theta(t)=\pi/2$ (or is even at minimum internuclear distance, as in the case when $\theta(0)=\pi/10$), then one observes that the molecule librates around $\pi/2$. However, since there is no coupling perpendicular to the field, the electronic forces along this component can only contribute to breaking the bond. In our model this is observed in $X_{2}^{(2)}(t)-X_{2}^{(1)}(t)$ quadratically increasing with time. We have not observed more than two re-crossings of the wave packet back and forth $\theta=\pi/2$ for any initial conditions, since the dissociation is relatively fast.

Figure 4 gives hints of the same physical mechanism, also providing information regarding the role of the electron. We show snapshots of the wave packet density at different times for the case $\theta(0)=\pi/10$. Initially the electron density is distributed evenly in both nuclei. As the bond relaxes, the electron moves with the proton at the right side ($x>0$), which is the expected result in the excited state in the presence of a strong field. Hence, as $\theta(t)$ crosses $\pi/2$ the electron density mostly shifts from one proton to the other, where the energy is negative. This is because in the excited state the density is larger at the nuclei sitting at positive values of the coupling.

IV Discussion

In principle, it is not difficult to determine the stereodynamics of photofragmentation by e.g. imaging techniquesChandler and Houston (1987); Eppink and Parker (1997).While very often one observes fragments perpendicular to the bond axis even in diatomic molecules, they are typically related to rotating molecules that give raise to isotropic signals, not to polarizability-induced anti-alignment effects. However, all the results in our simulations imply non-rotating molecules. The lack of experimental confirmation puts forth the question whether the results that we observe are an artifact of the approximations in our model. In the following discussion we provide compelling reasons to believe that the results are physical. As a first approximation, the interaction of the molecule with the field can be obtained from second order perturbation theory, which typically is also used to justify alignment. For the purpose of understanding the process it is enough to analyze the rotational dynamics, so we can consider a rigid diatomic molecule in a given electronic state. In nonresonant two-photon processes the effective interaction is given by Stapelfeldt and Seideman (2003); Seideman and Hamilton (2005)

[TABLE]

where $\alpha_{\parallel}$ and $\alpha_{\perp}$ are the parallel and perpendicular components of the polarizability with respect to the molecular axis. Clearly, the charge transfer resonance responsible for the strong binding (and large dipole) in H${}_{2}^{+}$, with transient dipoles increasing linearly with the internuclear distances and electronic state energies approaching asymptotically, dominates over any other dipole at large bond lengths, hence at large distances one only needs to account for $\mu_{12}\approx r/2$ along the molecular axis Chang et al. (2015) (the dipole between the ground and excited state, or the dipole between those states related by the charge resonance process). Then in the ground state, at large bond distance, takingAtkins (1983)

[TABLE]

where $\Delta(r)=(V_{2}(r)-V_{1}(r))$, is the detuning, or energy difference between the ground and first excited electronic states at the chosen internuclear distance, and we fixed the molecular axis in the $x$ direction. Then, $\alpha^{(2)}_{\parallel}\approx-\alpha^{(1)}_{\parallel}$ because of the changing sign in the denominator. Therefore, the same effect that leads to alignment in the ground state is responsible for the anti-alignment in the excited state. On the other hand, because of the large energy gap between $\Sigma$ and $\Pi$ states at large internuclear distances (or the 1s to 2p atomic energy difference), $\alpha^{(1)}_{\perp}\approx\alpha^{(2)}_{\perp}\approx 0$, so the perpendicular component of the polarizability cannot change the previous prediction.

In the excited state, the torque will be roughly given by $\tau^{(2)}=-dV_{\mathrm{ind}}/d\theta\approx\frac{1}{8}r^{2}E^{2}_{0}\sin(2\theta)/\Delta(r)$, so the aligned molecules are in a maximum of the potential, in an unstable configuration. That is, in the excited state the torque forces alignment of the molecular axis perpendicular to the field, where the field does not act and the molecule (in the excited state) dissociates. Obviously, this rough estimate using second order perturbation theory fails asymptotically when the electronic states are degenerate, but the calculation of the LIPs using perturbation theory for degenerate states shows the same effect. Note that in our results the calculation is exact (no perturbation theory is used) but the Hamiltonian is approximate (the electron is constrained to move in the molecular plane under a soft-core Coulomb potential).

V Summary and Conclusions

We solve the time-dependent Schrödinger equation for a $2$-dimensional soft-core Coulomb potential model of H${}_{2}^{+}$ in a grid under a strong field, using the Ehrenfest approach to account for the nuclear motion. The results for an ensemble of initial conditions reproduce qualitatively well the known alignment dynamics in the ground electronic state, and the creation of a huge oscillating dipole in the excited electronic state, when the molecules are initially aligned with the field. However, when the molecular axis is misaligned, we observe anti-alignment, that is, a torque that forces the molecular axis perpendicular to the field, where the dipole is zero, and the molecule inexorably dissociates.

The results of this work show that laser-induced vibrational trapping or stabilization in the excited state will be impossible in homonuclear cations, but we expect the presence of anti-alignment dynamics and its effects in other laser-induced processes to be very general, and prone to experimental measurement through e.g. imaging techniques, whenever a single transition dipole dominates at large internuclear distances between two states. On the other hand, if the photodissociation is the main goal of the dynamics, the anti-alignment process avoids laser-induced traps in the reaction coordinate. In addition, since at $\theta=\pi/2$ there is a LICI, the anti-alignment may help the wave packet to visit and hence to probe the dynamics in the vicinity of the LICI.

Acknowledgments

This work is supported by the Korean government through the Basic Science Research program (2017R1A2B1010215) and the National Creative Research Initiative Grant (NRF-2014R1A3A2030423), by the Spanish government through the MINECO Project No. CTQ2015-65033-P and FIS2016-77889-R, and by the Comunidad de Madrid through project Y2018/NMT-5028. FM acknowledges support from the “Severo Ochoa” Programme for Centres of Excellence in R&D (MINECO, Grant SEV-2016-0686) and the “María de Maeztu” Programme for Units of Excellence in R&D (MDM-2014-0377). BYC and IRS thank the Army Research Laboratory for the hospitality in a stay during which part of this work was created.