Magnetoelastic coupling in the cobalt adipate metal–organic framework from quasi-harmonic lattice dynamics

Abstract

Magnetic interactions in hybrid materials are poorly understood compared to those in purely inorganic materials. The high flexibility of many metal–organic systems introduces a strong temperature dependence of the magnetic exchange interactions owing to changes in the crystal structure. Here, we study the cobalt adipate system, for which anisotropic thermal expansion was recently shown to be a result of magnetoelastic coupling. The combination of density functional theory with quasi-harmonic lattice dynamics is shown to be a powerful tool for describing temperature dependent thermodynamic potentials that determine magnetic interactions. It is demonstrated that the effect of phonons can be sufficient to switch the preference for ferromagnetic versus antiferromagnetic ordering.

---

1 Introduction

Metal–organic frameworks (MOFs), consisting of small metal clusters connected by organic linkers, have received great interest due to the tunability inherent from the many possible combinations of metal and ligand.¹ Initial studies focused on porous structures because of their promising applications within gas storage² and catalysis.³ More recently attention has been drawn to properties typical of purely inorganic materials such as multiferroics^4–6 and extended magnetic ordering.⁷ These are often observed in dense frameworks which have shorter connection paths between the metal atoms than their porous counterparts.

The advantage of using MOFs as magnetic materials lies in their tunability. The magnetic coupling between two transition metal atoms is sensitive to the length and geometry of the super-exchange pathway, which can be varied by the rich choice of organic building blocks.⁸ A number of MOFs have already been synthesized with interesting magnetic properties such as low-dimensional magnetism,^8–10 spin-flop transitions,¹¹ magneto-caloric effects^12–14 and magnetoelastic coupling.^15,16 These effects are readily identified experimentally, however there is a need for a deeper theoretical understanding in order to direct further synthesis.

Here we use density functional theory (DFT) to study Co(adipate), Co(C₆H₈O₄), in order to understand the complex compromise between energy contributions from lattice geometry, magnetic interactions and lattice vibrations which lead to the observed magnetoelastic coupling. We show that the magnetic coupling as well as the phonon spectrum is dependent on the volume of the material. It is therefore necessary to consider the thermal expansion of the ferromagnetic (FM) and antiferromagnetic (AFM) states in the quasi-harmonic approximation to accurately predict the energy balance between them.

2 Computational details

Spin-polarised density functional theory calculations were performed using the Vienna ab initio simulation package (VASP)¹⁷ with PAW pseudopotentials. The PBEsol¹⁸ functional was used to describe the exchange–correlation energy. The magnetic unit cell (cf.Fig. 1), together with periodic boundary conditions, were used to model the infinite crystal, and the Brillouin zone was sampled with 1 × 2 × 2 k-points. For optimisation of the unit cell vectors at different volumes and for the calculation of lattice vibrations an energy cutoff of 700 eV was used, while a cutoff of 500 eV was used for the calculation of different magnetic configurations in the fixed AFM unit cell. As phonon calculations require tight convergence we furthermore set the EDIFF tag to 10⁻⁸, while 10⁻⁶ was used for all other calculations. All structures were relaxed to a maximum force of 0.01 eV per atom. Spin–orbit coupling was not considered, as the effect is generally weak for 3d-metals.

Fig. 1 Experimentally determined¹⁵ magnetic structure of cobalt adipate (a) in the ac- and (b) in the bc-planes. The cobalt atoms are coloured blue (magenta) for spin aligned parallel (antiparallel) to the a-axis. Carbon is grey, oxygen red and hydrogen white.

Lattice vibrations were calculated using Phonopy.¹⁹ The force constants were obtained by the finite-displacement method,²⁰ which in our case required force-evaluations for 114 symmetry-inequivalent displacements and fitting of the force vs. displacement to a harmonic function. The displacements were performed in the 1 × 1 × 1 magnetic unit cell due to the large size of this cell, and thus only the phonon dispersion at the Γ-point is considered.

In the harmonic approximation the bond lengths do not change with temperature, but the vibrational states calculated at the 0 K volume are populated according to Bose–Einstein statistics.

In the quasi-harmonic approximation^20,21 the effect of thermal expansion is taken into account by considering the anharmonic contribution arising from the volume-dependency of the phonon frequencies. The equilibrium volume (V) at a given temperature (T) and pressure (p) is the volume that minimises the Helmholtz free energy, F(T,V), to give the Gibbs free energy:

G(T,p) = min_V{F(T,V)}

= min_V{U(V) + F_phonon(T,V) + pV} (1)

where U(V) is the electronic energy calculated by DFT, and F_phonon(T,V) is the phonon energy at a given temperature and volume. We therefore calculate the phonon-frequencies in unit cells optimised at 8 different volumes and fit the Helmholtz free energy curves for a number of different temperatures to the Vinet equation of state²² in order to obtain G(T).

3 Results

The experimentally determined magnetic structure of Co(adipate)¹⁵ is shown in Fig. 1. The Co atoms interact antiferromagnetically via long-range super-exchange through the carboxylic acid groups of the adipate molecules in the b- and c-direction, while a weaker ferromagnetic coupling is found along the a-vector. Field-cooled magnetic susceptibility measurements show an onset of magnetic ordering at 15 K, which is accompanied by a change in the lattice parameters. The b- and c-vectors expand upon increasing temperatures, while negative thermal expansion is observed along the a-vector.

DFT was used to optimise the unit cell in the ferromagnetic (FM) and in the experimentally determined antiferromagnetic (AFM) state¹⁵ (cf. Computational section for details). The optimised lattice parameters are given in Table 1 along with the experimentally determined values at 10 K. The lattice parameters of the FM and AFM states are both close to the values determined by experiment, however small differences are predicted between the two states, resulting in the FM state having a slightly larger volume than the AFM state.

Table 1 DFT calculated lattice constants of the FM and AFM states, as well as the experimental lattice constants determined at 10 K

	DFT, AFM	DFT, FM	Exp.^15
a [Å]	15.99	15.93	16.147
b [Å]	9.67	9.69	9.586
c [Å]	9.42	9.46	9.248
β [°]	91.1	91.0	90.2
V [Å^3]	1455	1459	1432

---

The magnetic coupling can be described by the Ising-model Hamiltonian:

where J_ij is the exchange constant between neighbouring spins S_i and S_j. A positive value of J_ij indicates ferromagnetic coupling whereas a negative value indicates antiferromagnetic coupling.

The value of J_ij depends on the geometry of the exchange path. For the present system the exchange paths are identical along the b-vector, and can be described by a single exchange-constant, J_b. Similarly, the coupling is equivalent by symmetry in the c-direction and is described by J_c. In the a-direction two different couplings exist, but both are expected to be significantly smaller than J_b and J_c, and indistinguishable within the computational accuracy, and they are thus described by a single parameter, J_a.

In the simplest approximation we set J_a = 0, and assume J_b = J_c due to the similar distances between the Co-atoms along these two directions (4.83 Å along the b-direction and 4.72 Å along the c-direction). The exchange constant can then be found from the energies of the FM and the AFM states, E_FM and E_AFM, as:

where S_Co is the spin of one Co atom, N is the number of magnetic atoms in the unit cell and z is the number of interacting neighbours. Co is high spin in the tetrahedral configuration giving a value of , and each of the N = 8 ions interacts with z = 4 neighbours. Using the energies of the AFM and FM structures in their respective optimised unit cells results in an exchange constant of J = −0.08 meV. The uncertainty associated with comparing numbers calculated in different unit cells is estimated by calculating both E_AFM and E_FM in the same unit cell. For the optimised AFM unit cells this results in a value of J = −0.15 meV, while a similar calculation in the optimised FM unit cell gives J = −0.12 meV. We conclude that similar results are obtained by comparing numbers in slightly different unit cells, or by comparing numbers where either the AFM or the FM structure is not in its optimised unit cell. From mean-field theory it is possible to estimate the critical temperature from the exchange constant using the formula:which for the J-values given above results in critical temperatures in the range 10–18 K, in excellent agreement with the experimentally determined transition temperature of 15 K.¹⁵

We also attempted a more complicated analysis allowing different, nonzero values for all three exchange constants. The energy of 26 different magnetic configurations‡ were fitted to eqn (2) by a simple linear least squares method, giving values of J_a = 0.003 meV, J_b = 0.060 meV and J_c = −0.366 meV. This is not consistent with the experimental data, as it predicts ferromagnetic coupling along the b-vector, implying that the magnetic unit cell should be identical to the structural unit cell. The error is probably a result of the limited accuracy of DFT in describing weak magnetic interactions which involve very small energy differences. However, since the simple approximation in which J_b = J_c and J_a = 0 accurately reproduces the experimentally determined transition temperature the chosen method is still considered suitable for calculations involving the energy difference between those two states.

We now take a closer look at the relation between the magnetic interaction and the unit cell volume. To this end we optimise the FM and AFM structures at different volumes and calculate the exchange constant according to eqn (3) at each volume. The result, shown in Fig. 2, reveals that J increases with increasing volume, becoming positive (i.e. ferromagnetic) at a volume approximately 3% larger than the optimised volume of the AFM state (V_opt(AFM)). If, instead, the volume is decreased, e.g. as a result of an increased pressure, the AFM interaction becomes stronger and the transition temperature should increase.

Fig. 2 Exchange constant as a function of volume relative to the optimised volume of the AFM state, V_opt(AFM).

3.1 Phonon energies

The different unit cell dimensions and volumes of the FM and AFM structures mean that they will have different phonon energies. In the following we investigate how these differences affect the energy balance between the two states.

Fig. 3a shows the phonon partial density of states (PDOS) for the optimised AFM structure (cf. Computational section for details). The cobalt atoms are represented in the low-energy band below 500 cm⁻¹, while the range 500–1500 cm⁻¹ covers vibrations within the organic molecules involving the C-, O- and H-atoms. The isolated band at ca. 3000 cm⁻¹ contains the C–H stretch modes (see ESI† for .gif files of selected modes representing the three bands). The PDOS for the FM structure appears identical and is therefore not shown. Likewise, only subtle differences are found between PDOS obtained at different volumes. These are most clearly identified by plotting the zero point energy (ZPE) as a function of volume (Fig. 3b). The ZPE decreases with increasing volume for both the FM and the AFM states, and furthermore, even though the unit cells of the two states are not identical at a given volume, the ZPEs are nearly equal. The decrease in ZPE with increasing volume arises because the increasing bond lengths lead to slightly wider energy potentials, and thereby lower phonon energies.

Fig. 3 (a) Phonon partial density of states for the optimised AFM structure. Cobalt is shown in blue, oxygen red, carbon grey and hydrogen cyan. (b) Total zero point energy as a function of volume for the AFM (dashed blue line) and FM (full magenta line) states.

In the following we consider the difference in Gibbs free energy between the AFM and FM states as a function of temperature:

ΔG(T) = G_FM(T) − G_AFM(T) (5)

here G_FM and G_AFM are the Gibbs free energies of the system in the FM and AFM states, respectively, and a positive value of ΔG(T) thus indicates that the AFM state is more stable than the FM state.

Fig. 4 shows ΔG(T) with the phonon energies calculated in the harmonic approximation (red line) and in the quasi-harmonic approximation (black line). In the harmonic approximation the vibrational energy is calculated at the 0 K volume and the only effect of temperature is to excite the vibrational modes beyond the ground state. The free energy difference at 0 K, ΔG(0), is 6.8 meV which is smaller than the DFT internal energy difference of 12 meV. This is because the FM state has a larger volume than the AFM state, resulting in a lower ZPE. The energy difference becomes smaller with increasing temperatures, not as a result of an increasing difference in the vibrational energy, but rather as a result of the increasing contribution from vibrational entropy. The wider potentials associated with the larger volume of the FM state lead to larger vibrational entropy. As the entropy becomes more important with increasing temperatures the energy difference between the AFM and FM states decrease, however the AFM state remains the most stable over the shown temperature range.

Fig. 4 Free energy difference between the FM and AFM state as a function of temperature, with the phonons calculated in the harmonic approximation (dashed red line) and in the quasi-harmonic approximation (full black line).

In the quasi-harmonic approximation (Fig. 4, black line) thermal expansion is included by considering the volume-dependence of the phonon frequencies (cf. Computational section). This results in 0 K volumes that are approximately 3% larger than in the harmonic approximation, because the higher electronic energy of the larger unit cell is compensated by the softening of the phonon modes. The difference in electronic (DFT) energy at 0 K is 10 meV, but ΔG(0) is only 2.1 meV due to the difference in ZPE. As the temperature increases ΔG(T) becomes smaller and at 100 K the FM and AFM states are equal in energy.

3.2 Unit cell parameters

We now consider the effect that the magnetic transition predicted in the quasiharmonic approximation would have on the lattice parameters. Fig. 5 shows the lattice constants as a function of temperature for the AFM (blue) and FM (magenta) structures. The full lines indicate the lattice constants that would be observed assuming that the magnetic transition happens at 100 K, while the dashed lines indicate the lattice constants of the AFM (FM) state in their unstable region above (below) the transition temperature. The a-vector contracts when increasing the temperature above 100 K while the b- and c-vectors expand. Inspection of the dashed lines shows that the direction of change would be the same even if the transition happened at a different temperature. The same direction of change is observed experimentally when going from the ordered AFM to the disordered state. While we cannot model the disordered state by periodic boundary conditions DFT, we speculate that it has unit cell dimensions similar to the FM state, giving rise to the observed changes.

Fig. 5 Length of unit cell vectors for the AFM (blue) and FM (magenta) states as a function of temperature. The full lines show the unit cell parameters of the most stable magnetic state as predicted in the quasi-harmonic approximation, while the dashed lines indicate the parameters in the less stable state.

Our calculations demonstrate that the energy balance between the AFM and FM states is very sensitive, because the magnetoelastic coupling implies that a change in magnetic ordering is accompanied by a change in the vibrational energies as well as the vibrational entropy. In the present case we predict a transition from the AFM to the FM state at 100 K, however this is above the temperature where the transition to a disordered state is predicted to happen (10–18 K).

For the synthesis of a material where a change in magnetic state is possible, as well as for many other purposes, it would be desirable to have a higher critical temperature. This could be achieved by a stronger magnetic coupling, e.g. by use of a different metal such as Mn, which has 5 unpaired electrons in the high-spin state. Indeed, optimisation of the Co(adipate) structure with Co replaced by Mn leads to a predicted transition temperature of 96 K from eqn (4). This compound also exhibits magnetoelastic coupling, though in this case the a- and b-vectors change by less than 0.01 Å while the c-vector expands by 0.02 Å when going from the AFM to the FM state. Experimentally, however, the combination of Mn and adipate has been reported to result in two different crystal structures, both of which are different to the Co(adipate) structure.^23,24 The Mn₂(H₂O)(CO₂(CH₂)₄CO₂)₂ structure reported in ref. 24 consists of sixfold coordinated Mn atoms and has a low magnetic ordering temperature of 4.7 K, while half of the Mn atoms are sixfold coordinated and the other half are fivefold coordinated in the Mn₂(CO₂(CH₂)₄CO₂)(OH)₂ structure reported in ref. 23. Measurements of the magnetic ordering were not possible for the latter compound, but the related Mn₂(CO₂(CH₂)₂CO₂)(OH)₂ was reported to have an unusually high magnetic ordering temperature of 44 K, showing that this class of MOFs has potential as materials with strong magnetic coupling and interesting related properties.

Finally we note that magnetoelastic coupling has been observed in the formate perovskite [(CH₃)₂NH₂]Co(HCOO)₃.¹⁶ As in the cobalt adipate structure the metal atoms in this structure are connected through O–C–O bridges, suggesting this linkage as a motif of interest for further research into magnetoelastic materials.

4 Conclusion

In summary, we have studied the magnetic MOF cobalt adipate by first-principles lattice dynamics. The unit cells of the AFM and FM ordered states were found to be slightly different, confirming the experimental finding of magnetoelastic coupling. The exchange constant as well as the phonon energies were found to depend on the unit cell volume. This affects the energy balance between the two magnetic states, and if thermal expansion is taken into account the FM state becomes more stable than the AFM state above 100 K. The resulting change in the lattice constants is similar to the experimentally observed change upon loss of magnetic order.

Acknowledgements

We thank Jonathan M. Skelton for useful discussions. This work is supported by ERC programme grant no. 277757. The authors acknowledge computing support from the UK national supercomputing service (ARCHER), via membership of the UK Materials Chemistry Consortium which is funded by EPSRC (EP/L000202), and from the University of Bath computing services (Balena). PJS would like to thank the Glasstone Trust for funding via a fellowship.

References

1. G. Ferey, Chem. Soc. Rev., 2008, 37, 191–214.
2. S. Ma and H.-C. Zhou, Chem. Commun., 2010, 46, 44–53.
3. J. Lee, O. K. Farha, J. Roberts, K. A. Scheidt, S. T. Nguyen and J. T. Hupp, Chem. Soc. Rev., 2009, 38, 1450–1459.
4. P. Jain, V. Ramachandran, R. J. Clark, H. D. Zhou, B. H. Toby, N. S. Dalal, H. W. Kroto and A. K. Cheetham, J. Am. Chem. Soc., 2009, 131, 13625–13627.
5. A. Stroppa, P. Jain, P. Barone, M. Marsman, J. M. Perez-Mato, A. K. Cheetham, H. W. Kroto and S. Picozzi, Angew. Chem., Int. Ed., 2011, 50, 5847–5850.
6. G.-C. Xu, W. Zhang, X.-M. Ma, Y.-H. Chen, L. Zhang, H.-L. Cai, Z.-M. Wang, R.-G. Xiong and S. Gao, J. Am. Chem. Soc., 2011, 133, 14948–14951.
7. M. Kurmoo, Chem. Soc. Rev., 2009, 38, 1353–1379.
8. C.-B. Tian, C. He, Y.-H. Han, Q. Wei, Q.-P. Li, P. Lin and S.-W. Du, Inorg. Chem., 2015, 54, 2560–2571.
9. P. Canepa, Y. J. Chabal and T. Thonhauser, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 87, 094407.
10. P. J. Saines, J.-C. Tan, H. H.-M. Yeung, P. T. Barton and A. K. Cheetham, Dalton Trans., 2012, 41, 8585–8593.
11. P. J. Saines, P. T. Barton, P. Jain and A. K. Cheetham, CrystEngComm, 2012, 14, 2711–2720.
12. G. Lorusso, J. W. Sharples, E. Palacios, O. Roubeau, E. K. Brechin, R. Sessoli, A. Rossin, F. Tuna, E. J. L. McInnes, D. Collison and M. Evangelisti, Adv. Mater., 2013, 25, 4653–4656.
13. P. J. Saines, J. A. M. Paddison, P. M. M. Thygesen and M. G. Tucker, Mater. Horiz., 2015, 2, 528–535.
14. Y. Meng, Y.-C. Chen, Z.-M. Zhang, Z.-J. Lin and M.-L. Tong, Inorg. Chem., 2014, 53, 9052–9057.
15. P. J. Saines, P. T. Barton, M. Jura, K. S. Knight and A. K. Cheetham, Mater. Horiz., 2014, 1, 332.
16. R. I. Thomson, P. Jain, A. K. Cheetham and M. A. Carpenter, Phys. Rev. B: Condens. Matter Mater. Phys., 2012, 86, 214304.
17. G. Kresse and J. Hafner, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 558–561.
18. J. P. Perdew, A. Ruzsinszky, G. I. Csonka, O. A. Vydrov, G. E. Scuseria, L. A. Constantin, X. Zhou and K. Burke, Phys. Rev. Lett., 2008, 100, 136406.
19. A. Togo, F. Oba and I. Tanaka, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 78, 134106.
20. R. Stoffel, C. Wessel, M.-W. Lumey and R. Dronskowski, Angew. Chem., Int. Ed., 2010, 49, 5242–5266.
21. E. L. da Silva, J. M. Skelton, S. C. Parker and A. Walsh, Phys. Rev. B: Condens. Matter Mater. Phys., 2015, 91, 144107.
22. P. Vinet, J. R. Smith, J. Ferrante and J. H. Rose, Phys. Rev. B: Condens. Matter Mater. Phys., 1987, 35, 1945–1953.
23. P. J. Saines, P. Jain and A. K. Cheetham, Chem. Sci., 2011, 2, 1929–1939.
24. Y. Kim, E. Lee and D.-Y. Jung, Chem. Mater., 2001, 13, 2684–2690.

---

Footnotes

† Electronic supplementary information (ESI) available: .gif files of selected phonon modes. Optimised structures and output from phonon calculations can be retrieved from http://https://github.com/WMD-Bath/Phonons. See DOI: 10.1039/c5tc02633d

‡ The unit cell containing 8 Co atoms has 2⁸ = 256 possible magnetic states of which 28 are inequivalent by symmetry. Of these 26 were stable when optimised in the AFM unit cell, while the last two structures changed to different magnetic configurations during the relaxation.

---