Magnetic moment quenching in small Pd clusters in solution

Small palladium clusters in vacuum show pronounced magnetic moments. With the help of Born–Oppenheimer molecular dynamics simulations based on density functional theory, we investigate for the paradigmatic examples of the Pd13\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_{13}$$\end{document} and the Pd8\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$_8$$\end{document} cluster whether these magnetic moments prevail when the clusters are solvated. Our results show that the interaction with acetophenone quenches the magnetic moment. The reduction of the magnetic moment is a direct consequence of the electronic interaction between the Pd clusters and the solvent molecules, and not an indirect effect due to a different cluster geometry being stabilized by the solvation shell.

In this paper we study the magnetic moment of solvated palladium clusters. Our motivation for this is twofold and related to two of the above-mentioned possible practical uses of transition metal clusters. On the one hand, isolated small Pd clusters show large magnetic moments [20][21][22][23]25,27,28,[32][33][34][35]. For potential magnetic applications, it is of interest to know whether the magnetic moment prevails when the clusters are not in vacuum, but in interaction with surrounding media. This question is of interest also for a second reason. Pd nanoparticles and multimetallic clusters containing Pd have been demonstrated to be efficient in catalysis [29,30,36] and particularly noteworthy, excel in crosscoupling reactions [4,[37][38][39][40][41], as reviewed, e. g., in Refs. [42,43]. In several such experiments, the metal particles are embedded in some support, e. g., a metal organic framework, and the reaction takes place in solution. The question has been raised what really is the active species in such experiments [39,[44][45][46][47][48][49][50]: Is it the surface of the metal nanoparticle, or are very small clusters, e. g. of Pd, leaching out of the larger particle, i. e., go into solution, and form the catalytically active centers? For a review on this topic, see Ref. [51]. If the leached out particles retain a strong magnetic moment in solution, one might be able to experimentally answer the question whether leached particles contribute importantly to catalysis by detecting the magnetic moments and check for relations between magnetic and catalytic activity.
The phenomenon of leaching has been discussed for different solvents, e.g., water [50,[52][53][54] and organic solvents [44][45][46]. An accurate description of water with standard exchange-correlation (xc) approximations is a well-known challenge for density functional theory (DFT) [55]. Given the presence of d-electrons in the Pd clusters and the question of the reliability of the magnetic moment prediction of DFT in this context [56], studying the solvation of Pd clusters in water would be a task of considerable complexity. Our focus here is on acetophenone, which is the smallest aromatic ketone and an exemplary organic solvent whose derivatives have been used in studies of leaching and catalysis [44][45][46]. We investigated small Pd clusters in solution by computationally embedding them in an increasing number of solvent molecules and studying whether and how the magnetic moment changes as a consequence of the interaction. We find that acetophenone has a tendency to reduce the magnetic moment, i. e., a solvated Pd cluster shows a smaller magnetic moment than the same cluster in vacuum. We reach this conclusion based on extensive DFT-based Born-Oppenheimer Molecu-309 Page 2 of 13 Eur. Phys. J. D (2021) 75 :309 lar Dynamics (MD) simulations, the general setup of which is described in Sect. 2. In Sect. 3 we discuss how cluster structures and magnetic moments change when clusters are surrounded by an increasing number of acetophenone molecules. In Sect. 4 we discuss our results, also comparing to a few select findings for the solvent acetone that are reported in "Appendix B", and we summarize our results and draw final conclusions.

Theoretical and computational methods
The general concept of our study is to start from a Pd cluster with a nonzero spin magnetic moment in its vacuum geometry. We then add solvent molecules, find the low-energy structures that the system consisting of cluster and solvent molecules has, and check which magnetic moment is realized in these structures.
This procedure is straightforward in principle, but complicated in practice. First, we do not know a priori how many solvent molecules we need to take into account for the Pd cluster to be properly solvated. The computational effort escalates rapidly because, on the one hand, the computational costs of the DFT calculations increases with the number of electrons, and on the other hand, the number of low-energy isomers grows rapidly with the number of atoms. Therefore, we build up a solvation shell in a step-by-step manner. This means that we start with one acetophenone (C 8 H 8 O) molecule and advance to five, and more distant parts of the solvent are approximated via a conductor-like screening model (COSMO) as detailed below. Second, each system consisting of the cluster and several solvent molecules can be in many different geometries that are locally stable and of similar energy. Therefore, in order to identify low-energy structures, we first perform a Born-Oppenheimer DFT-MD simulation for some time and record the total energy. We then analyze this trajectory, pick out the lowest energy structures, and optimize those further by relaxing them.
Third, since the systems can assume different, stable spin configurations and since the spin configuration does not necessarily change spontaneously during a Born-Oppenheimer simulation, the just described steps have to be taken for each system size for all plausible spin configurations. The details of the procedure will be described in the following.
For our study we focused on the two small clusters Pd 13 and Pd 8 that show a distinct spin magnetic moment in vacuum [25,28]. They are furthermore good candidates because their vacuum structures have been established [10,22,23,25,27,28,32,57], and their isomers show canonical structural motifs like, e. g., the I h icosahedron for Pd 13 .
The energetic ordering of the isomers and spin configurations in clusters can generally depend on the xc functional [28,32,57]. For P 13 , the C s bilayer with 6 μ B shown in Fig. 1a and the C 3v bilayer with 6 μ B (b) is preferred by the Perdew-Burke-Ernzerhof [58][59][60] (PBE) and the PBE0 [61,62] hybrid functional, respec-  tively. For Pd 8 , the energetic ordering of the isomers has been discussed in detail in Ref. [63] and we thus limit ourselves to the Pd isomers relevant to our study. The D 2d bicapped octahedron shown in Fig. 2a is the previously reported global minimum for PBE [10,25,63,64]. It shows a 4 μ B ground state in our calculations. We find a state with 2 μ B just 0.02 eV higher in energy, i.e., quasi degenerate. With PBE0 the energetic difference between the two states increases to 0.10 eV. For the other two Pd 8 isomers (Fig. 2b, c) we also find 4 μ B as the preferred magnetic moment with both PBE and PBE0. In Tables 3 and 4 in "Appendix A", we list the energies for each of the Pd 13 and Pd 8 cluster isomers shown in Figs. 1 and 2, respectively, for both PBE and PBE0, for the low lying stable spin configurations.
For our calculations, we employed the Turbomole [65] program package. The geometry optimizations (relaxations) discussed in Sect. 3 used PBE in combination with the def2-TZVP [66,67] basis set and we used the effective core potential def2-ecp for Pd in combination with all def2 basis sets [66]. For the relaxation of each initial structure, we did several calculations with a fixed number of unpaired electrons to determine the overall lowest energy configurations. In other words, separate calculations were done for a range of spin magnetic moments μ. For, e. g., the ground-state relaxations of Pd 13 (C s ) and one acetophenone molecule we did calculations with 6 μ B , i.e., the moment of the bare metal cluster, and 0, 2, and 4 μ B . As mentioned in Sect. 3.1 we also did some calculation with 8 μ B as a crosscheck. In the DFT-MD simulations for the solvated clusters and the subsequent relaxation of selected low-energy structures, the studied magnetic moments could typically be restricted to the range of [0, 2, 4, 6]μ B due to the magnetic moment quenching that we discuss in detail below. Nevertheless, we conducted spot checks and also optimized structures with higher magnetic moments. This consistently confirmed that higher moments are energetically less favorable.
For simulating the structural dynamics of the solvated clusters we relied on constant temperature MD simulations with a Nosé-Hoover thermostat [68,69] set to a temperature of 298 K. Again the spin-state was fixed during each of these simulations. The Nosé-Hoover relaxation time τ = 560 a.u. ≈ 13 fs was seven times the length of the Born-Oppenheimer dynamics time step t = 80 a.u. ≈ 1.94 fs. Van-der-Waals' interactions were taken into account via the Grimme correction (DFT-D3) [70]. Solvent effects beyond the solvent molecules that we considered explicitly were modeled with the COSMO approach [71] using a relative permittivity r = 17.4 for acetophenone and r = 20.7 for acetone [72]. A total simulation time of ca. 6 ps was sufficient for obtaining reliable results about the dominant magnetic moments. For reasons of computational efficiency, we used PBE with def2-SVP [73] in the DFT-MD runs that served to identify low-energy geometries. The thus obtained low-energy structures were then further optimized with higher accuracy as described above.
We computed the binding energy E B of the composite M/S consisting of a metal cluster M and N solvent molecule(s) S, (1) as the difference of the molecular total energies E tot obtained from ground state DFT calculations, where E tot is the total energy for the system in its optimized geometry. Hence, E tot (M ) always refers to the structural isomer that provides the lowest energy with the functional employed, e. g., the Pd 13 (C s ) bilayer with 6 μ B for PBE. We indicate the xc functional used to calculate the total energies at the appropriate places throughout the article, unless it is clear that the results refer to the default choice PBE.
As part of our analysis of the interaction between clusters and solvent we checked for charge transfer between the cluster and the solvent molecules. For this we relied on the Bader charge density analysis [74] as implemented in the Bader code [75][76][77][78][79] to assign a charge to each of the atoms in the system. By comparing the charge of each atom in the composite system M/S to that within the separate subsystems M and S, we analyzed whether inter-molecular charge transfer occurs. We also analyzed the interaction based on the (generalized) Kohn-Sham density of states (DOS), obtained from the ground state eigenvalues. As a guide to the eye, we convoluted the eigenvalue spectrum with Gaussians of width 0.08 eV. For better inter-system comparison, we plotted the DOS on the same scale, [0:76] arb. unit, in all figures. In the calculations for the separate systems that are needed for both the Bader charge-and the DOS-analysis, we used the atomic coordinates that the subsystem atoms have in the relaxed, combined system to eliminate differences that would be due only to deviations from the vacuum structure.

Pd 13 with one acetophenone molecule
We first focus on the Pd 13 cluster because it is, as discussed above, a stable system with an established geometry. We start by looking at one acetophenone (0 μ B ) molecule capping a Pd 13 cluster. The C s bilayer structure, cf. Fig. 1a, with a magnetic moment of 6 μ B serves as the starting geometry in the following simulations.
Due to its geometry with a benzene ring and the CO and CH 3 groups, it is not a priori clear how the acetophenone molecule will orient itself with respect to the Pd 13 (C s ) bilayer cluster. Therefore, we performed several geometry optimizations starting from different arrangements of the acetophenone molecule as depicted in Fig. 3, and determined which combination of spatial orientation and magnetic moment is energetically favorable. In one case we aligned the acetophenone parallel to the cluster surface (Fig. 3a). In the other, acetophenone was set perpendicular with CO and CH 3 facing the Pd surface (Fig. 3b). For both cases, we started from several different positions of the acetophenone with respect to the cluster's surface, including inclined arrangements as shown in Fig. 3c and d. For all these arrangements, we ran geometry optimizations for different, fixed spin magnetic moments in the range of [0, 2, .., 8]μ B . We included the case of 8 μ B despite the fact that it would imply an increase in magnetic moment compared to the Pd 13 vacuum structure as a check. However, as expected the 8 μ B configuration is energetically clearly less favorable than the others.
The following observations were made in the relaxations: First, a spin configuration with 6 μ B is preferred for all orientations and in any case we do not observe  any significant distortions of the C s bilayer geometry through a single acetophenone. Second, in the case of an initially parallel orientation, our relaxed structures show that the benzene ring upholds the parallel orientation. In addition, the tilted parallel initial structures (cf. Fig. 3c) transition into a parallel arrangement (cf. Fig. 3a). For the parallel orientation, our binding energies range from −1.43 eV to −1.67 eV for a magnetic moment of 6 μ B . Third, relaxing the perpendicular initial situation leads to geometries in which CO stays closer to the C s bilayer surface than the CH 3 group, and the acetophenone remains perpendicular to the surface (similar to Fig. 3b). Inclined perpendicular initial structures (cf. Fig. 3d) generally transition into a clear parallel arrangement (cf. Fig. 3a). However, we also observed in some relaxations that initially perpendicular structures, both with and without inclination, flipped into a parallel geometry. In the perpendicular arrangement, the binding energies range from −0.27 eV to −0.78 eV (6 μ B ). One can thus conclude that the parallel orientation is clearly the preferred one. An explanation of this finding is provided by analyzing the spatial structure of the valence orbitals. For the parallel arrangement, their probability density delocalizes not only over the cluster, but to some extent also over the molecule, and this results in a lower energy. The magnetic moment of the combined system is the same as for the bare Pd 13 , but we occasionally observe stable configurations with lower magnetic moments. The first row of Table 1 summarizes the energetic ordering of the different spin configurations. In each case the table refers to the most stable structure that we found.
We also checked the robustness of the predicted spin configurations with respect to the xc approximation. To this end, we took the parallel initial structures that lead to the most stable configurations with PBE and repeated the optimization using PBE-D3 and PBE0-D3, respectively. The resulting binding energies and energy difference are listed in the second and third rows of Table 1. Reassuringly, both D3 corrected functionals also predict a 6 μ B spin magnetic ground state. The energetic ordering of the higher states changes to some extent, but this is not surprising given the small ener-  Table 5) the DOS of the combined system is dominated by Pd 13 (down to about −8 eV), whereas in the low energy region (below ca. 10 eV) it is dominated by C 8 H 8 O. The comparison of the DOS of the combined system to the sum of the two independent DOSs however reveals that there is also some electronic interaction, because the two are not identical. Calculating the DOS with PBE0 yields the same qualitative trend (cf. "Appendix A", Fig. 13).
The Bader charge density analysis for the parallel arrangement obtained with PBE-D3 shows a charge transfer of some 0.20 e − ≈ 0.2 e − from the Pd 13 cluster to the acetophenone molecule, where e − is the negative elementary charge. Performing the same analysis with PBE0 for the parallel arrangement obtained with PBE0-D3 yields the same qualitative trend with a charge transfer of some 0.26 e − ≈ 0.3 e − , i.e., the qualitative results are robust with respect to the influence of the xc approximation.

Pd 13 with three acetophenone molecules
We worked our way towards a more realistic description of solvation by combining the Pd 13 (C s ) cluster with three acetophenone molecules carrying out constant temperature DFT-MD simulations, following the procedure described in Sect. 2. We did separate DFT-MD runs for each of the magnetic moments μ ∈ [0, 2, 4, 6]μ B . As the start configuration we arranged two acetophenone molecules parallel and one perpendicular to the Pd surface. During the first ca. 1.2 ps the system restructures considerably as it moves from the initial geometry to a lower energy configuration: The one acetophenone that initially had a perpendicular orientation with respect to Pd 13 angles rapidly into a more stable parallel orientation, which further corroborates the results of the previous section. Then all three solvent molecules keep this type of orientation. This deforms the initial C s bilayer and for all of the following simulation time stabilizes a new Pd 13 structure with a higher spatial symmetry D 5h . The two basic structural motives are depicted by the insets in Fig. 5a.
In order to see which magnetic moment does now become the energetically preferred one, we record the total energy for each of the runs with different spin configurations. Because of the pronounced, rapid fluctuations in the total energy that are an inherent feature of MD simulations, the decisive trends are easier to see when one averages over some of those fluctuations. Therefore, Fig. 5 depicts the data in two different ways. First, we defined the average energyĒ(t) at a given time t = t i bȳ i. e., as the average over all of the previous discrete time steps i, including the current one. Here, 1 ≤ i ≤ N t and the total number of time steps was N t ≈ 3100. Averaging in this way over an increasing number of geometries allows to visualize the trend in the simulations, with fluctuations being totally averaged out. One thus obtains a direct depiction of the total energy trend for each spin configuration. (The more detailed information from a moving average is discussed below.) Figure 5a showsĒ ( Figure 5a shows that the DFT-MD run with zero magnetic moment leads to the lowest average energy. The 4 μ B structures are slightly higher in energy, and the 2 μ B and 6 μ B are noticeably higher. For a more detailed insight Fig. 5b depicts the total energy curve of the DFT-MD runs for each magnetic moment. In this plot we used a moving average E(t) of the total energy: At t i we averaged the total energy over the previous and subsequent Δ time steps: where Δ = 50, i − Δ ≥ 1 and i + Δ ≤ N t . Δ was chosen such that some of the smaller energy fluctuations were smoothed out, while preserving the main progression of the (unaveraged) total energy. Comparing E(t) for the different magnetic moments in Fig. 5b confirms that 0 μ B leads to the lowest energy.
In order to confirm yet further which magnetic moment is the energetically preferred one, we selected  Low-energy structures obtained from DFT-MD simulations, marked minMD = 1, . . . , 7 in Fig. 5b, served as the unrelaxed initial structures some of the low energy structures from each run. These are marked by the small circles numbered from 1 to 7 in Fig. 5b, and we refer to these numbers with the notation min MD in the following. We optimized (relaxed) these structures to the energetic minimum in separate calculations. Table 2 lists the binding energy E B (min MD ) for the relaxed structures. The relaxed structures show a 0 μ B or 2 μ B ground state. In some cases, the binding energies for the [0, 2, 4] μ B structures lie energetically close together, and one has to keep in mind that the employed xc approximations have a finite accuracy. However, the 6 μ B configurations are quite noticeably less favored. A geometrical analysis of the relaxed structures shows that the distance between the benzene rings of the solvent shell and the Pd surface is on average ca. 1.98Å.
Solvation thus clearly reduces the magnetic moment, and this raises the question of why that is. It could be a direct effect of the electronic interaction between the acetophenone molecules and the Pd 13 cluster, or it could be an indirect effect if the new D 5h geometry of Pd 13 that is triggered by the solvation is associated with a lower magnetic moment than the original C s structure. We checked for this by computing the magnetic moment of the D 5h structure in vacuum. First, we took the Pd 13 coordinates from the Pd 13 /(C 8 H 8 O) 3 arrangement and calculated the ground state for this bare Pd 13 structure without the acetophenone molecules. It shows a magnetic moment of 8 μ B . As a further check we also energetically optimized the extracted structure. The relaxed (bare) D 5h cluster also has an 8 μ B magnetic ground state, and it is 0.41 eV higher in energy than the C s bilayer structure with 6 μ B . Thus, the reduction of the magnetic moment is a direct consequence of the electronic interaction with the acetophenone.
The conclusion that there is a substantial interaction between the cluster and the solvent molecules is also confirmed by analyzing the DOS and the charge transfer. The DOS of Pd 13 /(C 8 H 8 O) 3 can neither quantitatively nor qualitatively be explained from the superposition of the DOSs of the separate subsystems Pd 13 and (C 8 H 8 O) 3 , as depicted in Fig. 6. The differences are very noticeable, e.g., in the peak structure close to the Fermi level. Compared to the scenario with a single acetophenone molecule, the differences are overall more pronounced. Also the charge transfer from the Pd 13 cluster to the solvent molecules is more pronounced

Pd 13 with five acetophenone molecules
By combining Pd 13 with five acetophenone molecules we reach a situation in which the metal cluster is surrounded by the solvent molecules. We used the Pd 13 (D 5h ) structure from the previous section, i. e., the one stabilized by three acetophenone molecules (cf. Fig. 1d) as the initial geometry for Pd 13 , and the five acetophenones were initially aligned parallel to the Pd surface. We then again ran separate DFT-MD simulations for each of the magnetic moments μ ∈ [0, 2, 4]μ B . Due to the results of the previous section we did not study higher magnetic moments. In the initial phase of the DFT-MD simulations we observed a restructuring of the central Pd 13 : The solvent shell stabilized an icosahedral I h cluster (cf. Fig. 7, inset). After this phase, the averaged total energyĒ(t) of the DFT-MD runs plotted in Fig. 7 shows a clear tendency towards a quenching of the spin magnetic moment to either 0 μ B or 2 μ B , with 4 μ B being noticeably higher in energy.
We extracted five of the lowest-energy structures from the DFT-MD runs, and their relaxations confirmed the trend: 0 μ B or 2 μ B are the preferred lowenergy states, while 4 μ B is clearly higher in energy by some 0.10 eV to 0.45 eV. The most stable structure,   We also tested again whether the reduction of the magnetic moment is a direct or an indirect effect by calculating the density and the magnetic moment for the bare Pd 13 in the geometry that it takes in the Pd 13 / (C 8 H 8 O) 5 arrangement. The bare extracted Pd 13 cluster recovers the 8 μ B that was also found in the optimized I h structure. This confirms the conclusion of the previous section that the quenching of the magnetic moment is a consequence of an electronic interaction between acetophenone and the metal cluster.
Again, the electronic interaction is also reflected in the charge density and the DOS of Pd 13 /(C 8 H 8 O) 5 : Firstly, the charge transfer from Pd 13 to the solvent shell continues to increase, to about 1.2 e − . The charge transfer occurs predominantly to the C atoms of the benzene-like ring of all five acetophenone molecules. Secondly, the DOS of the combined system clearly differs from that of the superposition of the stabilized Pd 13 cluster and the solvent shell (cf. Fig. 9), as seen, e. g., close to the Fermi level. Thus, the stabilized Pd 13 cluster interacts with the entire shell of five acetophenones.
An analysis of the orbitals of the combined Pd 13 / (C 8 H 8 O) 5 system with 0 μ B corroborates that there is   Fig. 10a for the highest-occupied molecular orbital (HOMO) of Pd 13 /(C 8 H 8 O) 5 with 0 μ B , computed with PBE. Because it is well known that semilocal exchange-correlation functionals tend to overestimate orbital delocalization, we repeated this calculation with the PBE0 functional. Figure 10b shows the HOMO from the PBE0 calculation, and the similarity to panel (a) is apparent. We have furthermore checked the 15 highest occupied orbitals and found that both functionals predict a qualitatively very similar delocalization for these. We can therefore conclude that the observed delocalization is not an artifact of the semilocal xc approximation.

Pd 8 with five acetophenone molecules
As the final step of our investigation, we check whether the reduction of the magnetic moment upon solvation is specific to Pd 13 , or also occurs for another system. To this end, we started from the C 2v bicapped octahedron Pd 8 shown in Fig. 2b with a 4  We then ran separate DFT-MD simulations of the composite system for each of the magnetic moments μ ∈ [0, 2, 4]μ B . In the initial phase of the DFT-MD simulations, we observed a restructuring of the central Pd 8 during which the solvent shell stabilized another bicapped octahedron isomer (cf. Fig. 11, inset) with a reduced symmetry C 1 . After this phase, the averaged total energyĒ(t) plotted in Fig. 11 shows a clear trend towards a quenching of the magnetic moment: Both 0 μ B and 2 μ B structures are clearly lower in energy than the 4 μ B ones. The relaxation of the two lowestenergy structures from the DFT-MD runs leads to a 0 μ B ground state and thus further confirms the reduction of the magnetic moment. The most stable arrangement, see Fig. 8c, has a binding energy of −5.98 eV with 0 μ B . 2 μ B is 0.17 eV higher in energy, and 4 μ B is clearly less favorable by 0.85 eV (cf. Table 1, line 5).
We again computed the magnetic moment for the metal cluster in the geometry that it takes in the Pd 8 / (C 8 H 8 O) 5 arrangement with 0 μ B . Without the acetophenone molecules, the magnetic moment goes up to 4 μ B , with a 2 μ B configuration just slightly higher (0.03 eV) in energy. We also again found a delocalization of the valence orbitals over the metal cluster and acetophenone molecules.
We also analyzed the charge density and the Kohn-Sham DOS in the same way as for the Pd 13 cluster surrounded by acetophenone. This revealed a charge transfer of about 1.0 e − from the Pd 8 cluster to the solvent shell. The DOS of the combined system Pd 8 /(C 8 H 8 O) 5 also differs noticeably from the DOS that results from summing the DOSs of the molecular components (cf. Fig 12). Thus, all the effects observed for Pd 13 are also found for Pd 8 , which is a strong hint that these findings are general.

Discussion and conclusion
We used DFT-Born-Oppenheimer MD simulations combined with further geometry optimizations to study how the magnetic moment of Pd clusters changes when those are solvated in the organic solvent acetophenone. Our study focused on Pd 8 and Pd 13 as the vacuum ground states of these systems are clearly established, have a nonzero magnetic moment, and as these clusters have isomers that show a range of paradigm geometries. Our simulations reveal a clear trend towards a reduction of the magnetic moment. For both clusters the calculations with the largest number of solvent molecules lead to ground states with a zero magnetic moment.
Analysis of the molecular orbitals shows that there is a pronounced interaction between the Pd particles and the acetophenone molecules, with the valence orbitals delocalizing over both the metal and the acetophenone. The interaction also manifests in a noticeable charge transfer from the metal cluster to the solvent molecules, and in clear signatures in the DOS. Putting these findings into perspective with earlier studies, we note that the presence of one benzyl alcohol and one benzylamine molecule also lead to a reduction of the magnetic moment of Pd clusters [18,19]. Futschek et al. reported that the quenching of the magnetic moment can also occur in the composite systems of Pd 13 with the ligands phosphine and thiol [25]. Interestingly, however, the quenching is less pronounced for phospines, as there the system retains a 4 μ B state even with high ligand coverage.
Our results show that detecting solvated Pd particles via their magnetic moment, and in this way contributing to answering the "leaching" question that was discussed in the introduction, is going to be difficult when the solvent is acetophenone. The situation may be different, though, for solvents that do not interact as strongly. Our study provides hints that the aromatic character of acetophenone plays an important role in the electronic interaction that leads to the magnetic quenching. In order to check this hypothesis, we computed as an outlook the interaction between Pd 13 and five acetone molecules. The results are reported in "Appendix B" and show that, indeed, the interaction and magnetic moment quenching are less pronounced for these smaller solvent molecules that do not feature an aromatic ring. Thus, avoiding the strong metalsolvent interaction triggered by the aromatic character is a close lying first step on the search for solvents that preserve the magnetic moment.

Author contributions
SH ran all calculations that are reported in the paper. SK supervised the work and SH and SK discussed the results regularly. LL introduced SH to technical aspects of using Turbomole. SH and SK wrote the manuscript, with suggestions and comments also from LL. All authors were involved in the scientific discussion of the results and in drawing the conclusions. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.

Appendix A: Electronic structure details
Tables 3 and 4 complement the discussion on the Pd 13 and Pd 8 cluster isomers of Sect. 2, respectively. Both tables list the total energy of the most stable arrangement (column 3) and the energy differences to configurations with other magnetic moments relative to the former (column 4 and the following) for PBE and PBE0 (column 2). Relaxing Pd 13 (D 5h ) with PBE0 did not lead to a stable configuration, and Pd 8 (D 2d ) with 0 μ B was not stable with PBE, therefore no numbers are reported for these cases.
As a further insight into the electronic structure we report how the binding between the Pd cluster and acetophenone affects the energetic position of the frontier orbitals. Table 5 lists the highest occupied molecular orbital (HOMO) eigenvalue and the size of the energy gap to the lowest unoccupied molecular orbital (LUMO) for the most stable arrangement that we found for each composite system M/S, and the same quantities for the different bare cluster isomers and solvent molecules.
For the most stable parallel arrangement of Pd 13 / C 8 H 8 O with 6 μ B obtained with PBE (line 1), the HOMO is at −4.46 eV and the gap to the LUMO is 0.10 eV. The comparison with the (bare) components Pd 13 (C s ) (line 6) and C 8 H 8 O (line 11) shows that both the HOMO and LUMO eigenvalue of the combined system are very similar to the values for the (bare) metal cluster. For the D3-corrected structures of Pd 13 /C 8 H 8 O and its components, the HOMO energy and gap deviate by only ±0.01 eV from the values given in Table 5 and are therefore not listed separately.
For the most stable arrangement of Pd 13 /(C 8 H 8 O) 3 , shown in Fig. 8a, the HOMO eigenvalue is −4.43 eV and the gap to the LUMO is 0.11 eV (cf. Table 5, Table 3 Total energy Etot (third column) for the PBE and PBE0 xc functional (xc, second column) for the lowest energy configuration for each of the Pd13 cluster isomers (first column) shown in Fig. 1  The numbers given refer to the most stable arrangement of each composite systems M /S (lines 1-5). For comparison, the same quantities are reported for the stabilized Pd clusters (lines 6-11) and the two solvent molecules (12)(13). Numbers are given in eV and refer to def2-TZVP/PBE line 2). The HOMO of the most stable arrangement of Pd 13 /(C 8 H 8 O) 5 (cf. Fig. 8b) is at −4.27 eV and the LUMO is 0.19 eV higher in energy (line 3). Thus, compared to the systems with one and three acetophenone molecules, the energy of the HOMO is slightly higher and the gap to the LUMO is slightly larger. The results show that the size of the acetophenone solvent shell has a rather small influence on the HOMO and LUMO. For completeness, we also included the most stable arrangement of Pd 8 /(C 8 H 8 O) 5 (cf. Fig. 8c) in Table 5. The comparison with the separate components shows that also here, the HOMO and LUMO energies are dominated by the Pd cluster (cf. lines 4, 10 and 11).
Finally, we also report results that allow to check the influence of the xc approximation. To this end we looked at the relaxed PBE0-D3 structures of Pd 13 / C 8 H 8 O, Pd 13 (C s ) and C 8 H 8 O. As expected, the hybrid PBE0 yields a lower HOMO energy and a larger HOMO-LUMO gap for each system: The HOMO of Pd 13 (C s )/C 8 H 8 O (Pd 13 ) is at −5.03 eV (−5.14 eV), and the LUMO is 1.69 eV (1.69 eV) higher. Otherwise, however, the qualitative picture closely corresponds to the one found with PBE. This conclusion is also supported by Fig. 13, which shows the PBE0 DOS computed within the relaxed PBE0-D3 structure for the parallel arrangement of Pd 13 /C 8 H 8 O, cf. Sect. 3.1.

Appendix B: Pd 13 with five acetone molecules
We here report results for a different solvent to provide an outlook on possible future work. A non-aromatic solvent makes for an interesting comparison, because our findings for acetophenone indicate that the aromatic ring plays an important role in the binding and interaction with the palladium clusters. The results reported here are to be understood as a first, qualitative test.
To keep the differences between the solvents small except for the missing aromatic ring, we remained in the family of ketones and selected acetone (C 3 H 6 O), which is the simplest of them and of similar polarity as acetophenone. We combined Pd 13 with five acetone molecules, which compares best with Pd 13 with three acetophenones in terms of the total system size, and with Pd 13 /(C 8 H 8 O) 5 in terms of the number of solvent molecules. We chose Pd 13 (I h ) as the initial geometry for the cluster and calculated optimized structures of Pd 13 / (C 3 H 6 O) 5 in analogy to the procedures described in the previous sections on acetophenone (cf. also Sect. 2).  The thus determined most stable arrangement, see Fig. 14, has a binding energy of −3.39 eV and prefers 6 μ B with PBE. An arrangement with 4 μ B is slightly less favorable by 0.04 eV, followed by 2 μ B by 0.08 eV, 0 μ B by 0.18 eV, and 8 μ B by 0.19 eV. The binding between Pd 13 and the acetone shell is thus smaller than in the acetophenone systems with similar system size or number of solvent molecules, respectively. For Pd 13 / (C 3 H 6 O) 5 , the quenching to 6 μ B -compared to the 8 μ B of the Pd 13 (I h ) cluster in vacuum -is also caused by an electronic interaction, as revealed by the analysis of the stabilized Pd 13 cluster, which we performed in analogy to the acetophenone systems. However, the quenching is considerably less pronounced than in the comparable acetophenone systems.
We then further characterized the electronic structure in analogy to the previous sections: For the most stable arrangement, the HOMO is at −3.86 eV and the HOMO-LUMO gap is 0.07 eV (cf. "Appendix A", Table 5, line 5). The comparison with the (bare) molecular subsystems Pd 13 (I h ) (line 7) and C 3 H 6 O (line 12) shows that both the HOMO and LUMO eigenvalues of the combined system are determined mostly by the (bare) metal cluster.
The binding between Pd 13 and (C 3 H 6 O) 5 is reflected in the corresponding DOS by qualitative and quantitative changes that one observes when comparing to the superposition of the molecular components, as depicted in Fig. 15. Overall, this suggests a direct electronic interaction between the cluster and the solvent, and in this respect the effects are similar to those found for acetophenone. However, the two solvents differ significantly with respect to the inter-molecular charge transfer: About 0.3 e − is transferred from Pd 13 to the acetone shell, which is a factor of about 2 and 4 less than for Pd 13 /(C 8 H 8 O) 3 and Pd 13 /(C 8 H 8 O) 5 , respectively. The reduced charge transfer, together with the reduced magnetic moment quenching indicate that nonaromatic solvents may better preserve the magnetic moment of solvated clusters.