Structure and spin of the low- and high-spin states of Fe²⁺(phen)₃ studied by x-ray scattering and emission spectroscopy

Electron transfer in transition metal complexes can convert photon energy into electric or chemical energy which has found applications in solar cells and photocatalysis.¹ In photosynthesis, the transfer of an electron from a metal to a ligand is responsible for the creation of a chemical potential that activates ATP, the energy carrier in biological systems.² In photocatalysis, the light driven splitting of water into dioxygen and dihydrogen, catalyzed by a photoactivated Ru(II) complex, is an example of the strive toward cleaner and renewable energy.³ 

In transition metal complexes with open 3d, 4d, or 5d shells, a d electron is transferred upon photon absorption to a ligand orbital, and this state is called the metal-to-ligand charge transfer (MLCT) state. The process is reversible and the electron returns to the d orbital of the metal. In 4d and 5d complexes, the LMCT lifetime is long, from nanoseconds to microseconds. In 3d complexes it is of the order of sub-picoseconds in most cases. In 4d and 5d complexes, the outer electron density is predominantly 4d or 5d in character and the ligand bonding is stronger than for 3d complexes. The ligand field splitting of 4d or 5d orbitals is therefore greater than for 3d complexes. In fact, the lowest excited state for 4d and 5d complexes is the LMCT state whereas for 3d complexes it is a metal centered (MC) ligand field state. In the decay of the MLCT state in 3d complexes, the electron returns to the ground state (GS) via the MC state which greatly shortens the LMCT lifetime compared with the equivalent 4d and 5d systems.⁴ In the 4d and 5d complexes Ru²⁺(bpy)₃ and Os²⁺(bpy)₃, the MLCT lifetimes are 900 (Ref. 5) and 25 ns (Ref. 6), respectively, whereas it is 0.96 ps for Fe²⁺(bpy)₃ .⁷ The short MLCT lifetime of 3d complexes makes them unsuitable for light conversion since the charge separated state does not have time to interact with the environment. Since 4d and 5d metals are rare and expensive, it is important to synthetize ligands for 3d systems with the strongest possible ligand field. The ligand field strength depends on the electronegativity of the ligand atoms in the metal–ligand bonds: O and halogens produce weak fields, N intermediate fields, and C and P strong fields. In the carbene complex Fe³⁺(btz)₃(PF₆)₂, the strong Fe–C σ bonds generate a high electron density on Fe³⁺ and the MC states are inaccessible in the decay process. The MLCT state decays directly to the GS by photon emission with a long 528 ps lifetime.⁸ Another example is the Fe³⁺(phtmeimb)₂ complex, which has a record MLCT lifetime of 2 ns.⁹ 

In this work we present a time resolved study of the change in structure and spin of the HS state of Fe²⁺(phen)₃ (phen = 1,10-phenanthroline) in water using wide-angle x-ray scattering (TRWAXS) and x-ray emission spectroscopy (TRXES). The aim is to measure the change in the average Fe–N bond length to the ligands, the relative ligand angles, the spin states of Fe as well as the cage structure and hydrodynamics of the solvent. The work was performed by laser/x-ray pump probe using 100 ps x-ray pulses from the European Synchrotron Radiation Facility (ESRF).^10,11 The experimental details are given in Sec. II and in the supplementary material.³⁵ 

Fe²⁺(phen)₃ has the stoichiometry Fe²⁺(N₂C₁₂H₈)₃ with 67 atoms per molecule. Phenanthroline, a medium strong ligand, has two Fe–N bonds per ligand in Fe²⁺(phen)₃. The LS and HS structures were calculated by density functional theory (DFT) using the B3LYP function and the 6-311++G(3df,3pd) basis set. The effect of the solvent was included using the polarizable continuum model (CPCM) with Fe²⁺(phen)₃ placed in a cavity in the solvent field. The LS structure was assumed to have D₃ symmetry in the calculation. The LS calculation predicts six equal Fe–N bond lengths of 1.99 Å and relative ligand angles of 88.8°. The HS state was assumed to have C₂ symmetry. The Fe–N bonds and the relative ligand angles were calculated to be 2.18 Å and 94.2° on average. The LS structure is shown in Fig. 1(a). The atomic positions and structural parameters are listed in Appendix 1 in the supplementary material.³⁵ 

The 3d orbitals in the ligand field states t_2g and e_g are shown in Fig. 1(b) and the LS and HS configurations are shown in Fig. 1(c).

In the ground state t_2g, the lopes of the 3d_xy, 3d_xz, and 3d_yz wavefunctions are pointing along the diagonals in the (x, y), (y, z), and (z, x) planes towards the midpoints of the N-N lines of the ligands where the repulsive energy is lowest. In e_g the d_x2-y2 and d_z2 lopes are pointing along the x, y, and z axes in the direction of the N atoms. The repulsive energy is higher and the Fe–N bonds are longer.

The electron configuration of Fe²⁺ is [Ar]3d⁶. The 3d orbitals are tenfold degenerate and can accommodate five up-spins and five down-spins. The six 3d electrons in Fe²⁺are distributed in t_2g and e_g according to Hund's rule: maximize the total spin at the lowest possible total energy. If the ligand field splitting is greater than the spin-pairing energy in an orbital, the case for Fe²⁺(phen)₃, the three orbitals in t_2g are filled with three up- and three down-spins. The total spin is S = 0. For the next 7–10 electrons in the following 3d metals, e_g is filled by maximizing the spin. For ten 3d electrons, the total spin S is  0. If the ligand field splitting is smaller than the spin-pairing energy, the first five electrons have parallel spins in t_2g and e_g. The next five electrons are then put in t_2g and e_g with antiparallel spins.

As mentioned, phenanthroline is a medium strong ligand and the GS of Fe²⁺ (phen)₃ has six 3d electrons in t_2g with antiparallel spins and S = 0. When the e_g orbitals are populated in HS at the end of the photocycle, e_g receives the recombining electron which triggers the transfer of one t_2g electron to e_g with a change in spin. The two t_2g electrons have parallel spins and the total spin S is  2.

The photocycle of Fe²⁺(phen)₃ in acetonitrile was studied by ultrafast optical spectroscopy by Tribollet et al.¹² and Gallé et al.,¹³ by time resolved EXAFS in aqueous solution by Nozawa et al.,¹⁴ and by XANES by Wang et al.¹⁵ The Fe–N bond lengths derived from the EXAFS work are 1.98(1) and 2.15(1) Å for the LS and HS states, respectively. The HS lifetime was reported to be 690 and 685 ps in the EXAFS and XANES works, respectively. For comparison, the HS lifetimes of the sister complexes Fe²⁺(terpy)₂ and Fe²⁺(bpy)₂ in aqueous solutions are 2.61 ns (Ref. 34) and 0.665 ns (Ref. 7), respectively, and the Fe²⁺(btbip)₂ complex has a HS lifetime in acetonenitrile of 0.23 ns (Ref. 16).

The Fe²⁺(phen)₃ wavefunctions, transition rates, and energy curves were calculated in acetonitrile by Sousa et al.¹⁷ and the energy curves and spin states are summarized in Fig. 2. The photocycle proceeds as follows. When Fe²⁺(phen)₃ is excited from the GS ¹A₁, a 3d electron in t_2g moves to one of the ligands to form the ¹MLCT state. The ¹MLCT transits in 36 fs to the triplet ³MLCT that in turn decays to the MC triplet ³T₁ in 150 fs. ³T₁ decays in turn to ⁵T₂ in 102 fs by intersystem crossing (ISC). The ⁵T₂ state is reached in ∼300 fs. ⁵T₂ is initially in an excited vibrational state that relaxes to the bottom of the potential in a few picoseconds via vibrational cooling by the solvent. The final step ⁵T₂ → ¹A₁ proceeds via ISC in 1.1 ns in acetonitrile. Note that ³T₁ and ⁵T₂ have one and two antibonding electrons in e_g, respectively, which explains the progressive elongation of the Fe–N bonds.

The HS(S = 2) state was first observed by Tribollet et al.¹² Their Fe²⁺(phen)₃/acetonitrile solution was excited by 30 ps pulses at 540 nm and probed by 30 ps pulses of white light. They observed a transient absorption feature around 340 nm with a 1.1 ns lifetime that they assigned to the HS state ⁵T₂. The absorbing feature was assigned to the ⁵T₂ → ⁵MLCT excitation.

TRWAXS is used in this work to monitor the structural changes in the LS to HS transition, i.e. the change in Fe–N bond lengths, the relative ligand angles, the effective radius of the solvent cage, and the change in the solvent temperature. As 1D scattering data cannot be inverted to electron density directly, the measured scattering curves were compared with the predictions from combined DFT and MD simulations. The change in the spin state of Fe²⁺ was probed by TRXES from the K_β emission line. The experiments were performed on beamline ID09 at the European Synchrotron Radiation Facility (ESRF) and the time resolution was obtained by pump and probe. A 10 mM aqueous solution of Fe²⁺(phen)₃ was circulated in a liquid jet, excited by 1.2 ps pulses at 400 nm and probed by 100 ps x-ray pulses from an undulator. In the TRWAXS study, 18 keV x-ray pulses were produced by a 2.1% bandwidth (BW) Ru multilayer (ML) to maximize the intensity while keeping the beam quasi monochromatic. The intensity of the ML beam was 4 × 10⁸ ph/pulse and the pump & probe pulses were repeated at 986.3 Hz. In the TRXES experiment, the 5% BW pink beam from the fundamental of the U17 undulator at 15 keV was used to generate 1s core holes on the Fe K-edge. The intensity of the pink beam was 2 × 10⁹ ph/pulse. With the 100 ps time resolution from the x-ray pulse length, the HS state is the only structure resolved in the present experiment. The setup and the experimental parameters are described in more detail in the supplementary material.³⁵ 

The scattering functions S(q) for the LS and HS structures are shown in Fig. 3(a) and qΔS(q) are shown in Fig. 3(b). q multiplication is used to amplify the high-q spectrum, which is particularly sensitive to intra molecular structural changes. Note that qΔS(q) has 12 nodal points in the range 0 < q < 10.7 Å⁻¹, which is a very stringent constraint for the comparison of calculated and measured scattering curves.

The measured curves ΔS(q, t) have two principal contributions: one from the change in the solute/cage structure and one from the change in the solvent temperature from the relaxation of the photoproducts. The measured qΔS(q, t) curves for 100 ps and 2 ns are shown in Fig. 4. The blue curve in Fig. 4(a) is the total change in the solute/cage/solvent scattering (HS-LS). The solvent part is produced by the adiabatic temperature rise in the solvent after 100 ps. The orange curve is the solvent contribution that was measured in the separate dye/water experiment described in the supplementary material.³⁵ The difference between the curves is the intrinsic change in the solute and solute cage structure. The negative part between 0.5 and 1.0 Å⁻¹ in the blue curve is the signature of the Fe-N bond extension in HS. The rise in the solvent temperature after 100 ps is produced by the ultrafast cascade ¹MLCT → ³MLCT → ³T₂ → ⁵T₂ shown in Fig. 2. For the 2 ns curves in Fig. 4(b), the total solute/cage/solvent curve in blue coincides with the dye/water curve: the HS state has returned to the GS and the residual signal is from the adiabatic temperature rise of the solvent. Note that at fast timescales t < 10 ns, the solvent temperature rises monotonically and that the rise accumulates all relaxation processes in the photocycle, including the ultrafast ones that are not resolved structurally with 100 ps time resolution.

The normalization of the intensity is done by scaling the non-excited total scattering to 1.0 at q₀ = 6.7 Å⁻¹ and using that calibration for the difference scattering. The total scattering from a 10 mM “unit-cell” solution of Fe²⁺(phen)₃ in water, elastic plus Compton scattering, is 13.3 (e.u.)²/unit.^20,21

The MD simulations of the solvent structure around LS and HS states were performed by MOLDY employing the Lennard–Jones potentials from the TIP4P/2005 model of water.²⁵ In the simulation, the LS and HS structures of Fe²⁺(phen)3 were placed in a (49.6 Å)³ box with 4096 water molecules. The atom–atom distribution function g_Fe,Ow(r), where O_W is oxygen in water, is a measure of the effective cage radius. The g_Fe,Ow(r) and g_N,Ow(r) functions are shown in Fig. 5.

The first peak in the LS function g_Fe,Ow(r) at 5.5 Å shifts to ∼4.5 Å in HS. The latter has three sub-peaks separated by about 0.4 Å. The splitting is caused by the distortion of the HS structure, in particular the slightly different relative ligand angles in HS that vary from 90.4° to 96.0°. When the Fe²⁺(phen)₃ structure is rendered with 3D van der Waals spheres, the increased solvent accessibility to Fe in HS is evident, see Fig. 6. This accounts for the counterintuitive fact that water moves closer to Fe in HS in spite of the Fe–N elongation.

A splitting of the g_Fe,Ow(r) function in HS was also observed in MD simulations of Fe²⁺(bpy)₃ in water by the Hauser group²⁶ although the effect is less pronounced.

To optimize the fit of the experimental data to the DFT/MD model, the amplitude of ΔS_cage(q) is reduced to 33% of the MD calculated value. MD tends to overestimate the cage amplitude which is due to the fact that MD simulations are using rigid solute structures that allows the solvent to get closer to the solute.

The amplitude of the total scattering curves ΔS(q,t) decreases exponentially with a 725 ps time constant. The shape of the spectra vs. time does not change which is consistent with the fact the HS → LS transition is ultrafast for single molecules. By contrast, the amplitude decreases vs time from the decay of the population of the HS state.

Note that the cage signal is greatest between 0-2.5 Å^−1. and that the positive part signals the contraction of the cage. In the low q range 0.3 – 0.8 Å⁻¹, a small shift dr in the atom–atom distribution leads to a proportional change in the amplitude of ΔS(q) as shown in the supplementary material.³⁵ In that approximation, the cage radius contracts by 0.33 Å in response to the metal–ligand elongation in the HS state.

The cage dynamics for Fe²⁺(bpy)₃ in water was studied by TRWAXS by Khakhulin et al.²⁷ and Haldrup et al.²⁸ In Khakhulin's work, the solute and cage were calculated by ab initio MD (AIMD), which is a more realistic account of the solute/cage structure in water. In this complex the Fe–N bond extends from 1.99 Å (LS) to 2.18 Å(HS) and the AIMD calculation showed that the first peak in g(Fe, O_w)moves from 5.29 Å in LS to 5.51 Å in HS corresponding to a 0.22 Å expansion of the cage. The difference in the cage dynamics in Fe²⁺(phen)₃ and Fe²⁺(bpy)₃ is ascribed to the difference in the geometry of the ligands and the different structural distortions in the HS states.

As the temperature term dominates the dS(q,t) curves it is essential to have accurate solvent differentials to isolate the solute/cage contribution from the thermal background. This was done in a separate water/dye experiment performed with the same experimental conditions as the main experiment. The dye experiment is described in the supplementary material.³⁵ 

The slow rise ΔT₁ = 34 mK is consistent with the expected rise from a 2.4 mM excited solution at time zero and that the ⁵T₂ → GS decay deposits 0.62 eV of energy to the solvent. The excited state concentration derived from ΔT₀, which was expected to deposit 2.48 eV (3.10–0.62) of energy per molecule, is 5.1 mM. If the photoexcitation was a single-excitation event in this experiment, the two concentrations should be equal. We will now examine the discrepancy in the temperature rise since this effect has been observed in other complexes as well.

The excess solvent temperature rise might be explained by the fact that the 400 nm 1.2 ps pump pulse re-excites nacent ⁵T₂^*(HS) molecules to the ⁵MLCT band during the vibrational relaxation toward the bottom of the ⁵T₂ potential. Note that the formation of the ⁵T₂ state was indeed monitored by the ⁵T₂ → ⁵MLCT transition at 340 nm (3.65 eV) in acetonitrile by Tribollet et al.¹² Our 400 nm (3.10 eV) pump pulse might therefore excite ⁵T₂^* into the ⁵MLCT band in the early stage of the vibrational relaxation process. Note that the excitation-during-relaxation technique was used by C. B. Harris et al. in their classical work that probed the vibrational relaxation of I₂^* in CCl₄ and in other solvents.²⁹ The ⁵T₂* → ⁵MLCT transition is spin allowed and the transient absorption should increase for a 400 nm (3.10 eV) pump pulse once the molecule has relaxed to match the 3.10 eV excitation window. The ⁵MLCT state might then transit to ¹MLCT by ISC whereby the photocycle resumes. The re-excitation pathway is illustrated in Fig. 9(b).

To examine this hypothesis, an optical–optical experiment was performed in Mathias Bargheer's laboratory at the University of Potsdam. A 2 mM water solution of Fe²⁺(phen)₃ was excited by 400 nm 150 fs pump pulses and probed with 150 fs pulses at 400 nm and 500 nm. The results are shown in Fig. 9(a).

When the solution is probed at 500 nm near the 512 nm absorption peak, the absorption decreases with time from the depletion of the GS. For the 400 nm probe however, a Gaussian-like peak with a FWHM of 0.36 ps appears 0.19 ps after the pump pulse. The peak is tentatively assigned to the transient absorption from the ⁵T₂* → ⁵MLCT excitation with ⁵MLCT returning to ¹MLCT by ISC. A second photocycle turn is thereby initiated. Double-excited molecules dissipate 2 × 3.1 eV to the solvent, triple-excited 3 × 3.1 eV and so on.

For the 400 nm probe, the amplitude of the depletion hole and the creation peak is −0.064 and +0.098, respectively, which means that the cross section σ(⁵T₂* → ⁵MLCT) is 53% greater than σ(¹A₁ → ¹MLCT). That shows the importance of the proposed re-excitation pathway.

While WAXS probes atom–atom distances, intra- and extra-molecular alike, WAXS is insensitive to the spin of the molecules. In contrast, XES is atom-specific and a powerful probe of the spin and valence states of transition metals. As shown below the XES signals from metal complexes in solution are feeble but essentially background free. XES is used here to show that the spin of Fe²⁺ is 0 and 2 in the LS and HS states respectively.

Upon the absorption of an x-ray photon at the K-edge of Fe²⁺, a 1s electron is ejected and leaving behind a 1s core hole. The core hole is filled by an outer shell electron accompanied by the emission of either an x-ray photon or an Auger electron. For the K-edge of Fe, the branching ratio is 34% x-ray emission and 66% electron emission.³⁰ For the K_α and K_β lines, the emission is produced by the transitions $2 p → 1 s$ and $3 p → 1 s$⁠, respectively. The spin–orbit coupling of the 2p state splits K_α into two sub-lines K_α2 and K_α1 from the transitions 2p_1/2 (2e-) → 1s and 2p_3/2 (4e-) → 1s with the parenthesis showing the electron occupancy in that state. The K_α2 and K_α1 lines peak at 6390.8 and 6403.8 eV, respectively. The energy splitting of the K_β emission from 3p_1/2 (2e-) → 1s (K_β′) and 3p_3/2 (4e-) → 1s (K_β1,3) is dominated by the exchange interaction between the 3p core hole and the 3d electrons. The spin–orbit interaction for Kβ is much weaker than for K_α, about 0.7 eV (Ref. 31) due to the greater shielding of the nucleus. The central energy of K_β is 7058 eV for Fe.

XES is insensitive to the solvent hydrodynamics and the technique is ideal for monitoring population kinetics. One nice feature of XES compared with other spectroscopies is that the core holes can be produced by the pink beam which is ∼200 times stronger than the beam from a Si monochromator. The high beam intensity of 2 × 10⁹ ph/pulse from the undulator was essential in the present experiment due the low pump–probe frequency of 986.3 Hz.

The spin sensitivity of the K_β emission is due to the strong exchange interaction between the 3p core holes and the 3d electrons. The different spin states of Fe, singlet, doublet, …, quintet, have a unique K_β spectrum which is essentially independent of the ligand and the solvent.³² The valence-to-core (VtC) transitions 3d → 1s probe the ligand-modified 3d levels and VtC can therefore provide valuable information about the metal-to-ligand bonding.³³ However VtC is two-orders of magnitude weaker than K_β and too weak to be used here due to the low pump–probe frequency (986.3 Hz).

The K_β spectra ON and OFFand their HS-LS difference are shown in Fig. 10.

The shape of the difference spectrum HS(100 ps)−LS has the “pure” shape of the LS (S = 0) → HS (S = 2) transition since the non-excited molecules cancel out in the difference spectrum. The amplitude is proportional to the degree of excitation. The ON spectrum is shifted +0.2 eV from the OFF peak at 7058.3 eV. The appearance of the K_β^′ peak at 7044.6 eV is due to the increase in exchange splitting of the 3p level in the magnetic HS state. The VtC peak at 7110 eV is also detected; it is 100 times weaker than K_β.

The K′_β and K_β1,3 peaks in the HS spectrum are at 7045 and 7061 eV, respectively, corresponding to an exchange/spin–orbit splitting of 16 eV. The spectra in Fig. 11 were fitted to theoretical LS and HS spectra calculated by the CRISPY software. To match the experimental data, the calculated spectra were shifted +53 eV and broadened by 1.6 eV Lorentzians. The XES spectra from Fe²⁺(phen)₃ in water resemble closely those for Fe²⁺(terpy)₂ in water as reported by G. Vanko et al.³⁴ 

The photon-induced LS to HS transition of Fe²⁺(phen)₃ in water was probed by TRWAXS and TRXES using synchrotron radiation with 100 ps resolution. The TRWAXS data, measured over the wide q range 0.5–10.7 Å⁻¹, allowed to probe multiple length scales, from atom–atom distances inside the complex, the solvation structure around it, and the hydrodynamics of the solvent.

The difference scattering ΔS(q, t) is dominated by the change in the solvent temperature from the relaxation of the ¹MLCT (< 1 ps) and HS (725 ps) states. After thermal correction of the scattering data, the HS and LS structures were compared with DFT/MD simulations that provided trial functions for the LS and HS states and their cage structures. The high-q range 4.0–10.7 Å⁻¹ is in excellent agreement with the calculated Debye scattering from DFT with a 0.19 Å Fe–N bond elongation in the LS → HS transition. The elongation is caused by the antibonding character in the HS state (t_2g)⁴(e_g)² that has two 3d electrons in the antibonding orbitals d_x2-y2 and d_z2.

The solvent cage around the LS and HS states were simulated by MD and the S_cage(q) term was calculated from the g_i,j(r) functions. MD predicts a 1.0 Å contraction of the first solvation shell and a positive cage amplitude in the low-q range 0–0.8 Å⁻¹. The experiment shows that the cage amplitude is overestimated by MD and by reducing it to the 33% level, the fit to the data greatly improves. The amplitude reduction in q space corresponds to a 0.33 Å contraction of the cage radius in the HS state.

The time dependence of the solvent temperature was derived from the thermal solvent differentials from a separate dye/water experiment. The temperature rise from the slow 725 ps HS→LS relaxation is well explained by the energy difference of 0.61 eV between the HS and LS states and the initial HS concentration derived from the Debye fitting. In contrast, the temperature rise from the ultrafast (<1 ps) ¹MLCT → HS transition is 2.1× greater than expected from the DFT calculated energy difference of 2.49 eV (3.10–0.61) and the HS concentration. The anomaly is caused by the re-excitation of nascent HS molecules by the 1.2 ps pump pulse at 400 nm. . Ultrafast optical/optical absorption measurements showed the appearance of a new absorption feature after 0.19 ps that is tentatively assigned to ⁵T₂ → ⁵MLCT. ISC between ⁵MLCT and ¹MLCT then provides a link that restarts the photocycle which explains the excess temperature rise of the solvent.

Finally, the spin states of LS and HS were probed by XES with the core-holes produced by a pink undulator beam. When the 100 ps difference spectrum is compared with reference spectra for putative spin changes, it is assigned to the S = 0 to S = 2 transition. In the HS state, the K_β1,3 peak from the 3p_3/2 → 1s transition has a weak side-peak K_β′ from the transition 3p_1/2 → 1s which is −16 eV from K_β1,3 . The increase in the K_β′ intensity is due to the greater exchange interaction in the magnetic HS state.

We would like to thank Professor Carmen Sousa from the University of Barcelona most sincerely for explaining her theoretical work for us in great detail. We are also indebted to Dmitry Khakhulin and Kristoffer Haldrup for sharing their insight and experience with us. We also wish to thank Denis Leshchev, Norman Kretzschmar and Martin Nors Pedersen for help during the early stages of the project. We would also like to thank Laurent Eybert and Wolfgang Reichenbach for their work on the experimental setups on beamline ID09 at ESRF. The use of the femtosecond transient absorption spectrometer in Mathias Bargheer's lab at the University of Potsdam is also gratefully acknowledged. Finally, we are grateful for the EU Grant X-PROBE (Grant 637295) that funded Victoria Kabanova's Ph.D. studentship.

The authors have no conflicts to disclose.

Victoria Kabanova: Investigation (lead). Mathias Sander: Investigation (supporting). Matteo Levantino: Investigation (supporting). Qingyu Kong: Formal analysis (equal). Sophie Canton: Investigation (supporting). Marius Retegan: Formal analysis (supporting). Marco Cammarata: Validation (equal). Philipp Lenzen: Validation (supporting). Latévi Max Daku Lawson: Methodology (supporting). Michael Wulff: Investigation (equal); Validation (equal); Writing – original draft (equal).

The data that support the findings of this study are available within the supplementary material.