H3O+ tetrahedron induction in large negative linear compressibility

Despite the rarity, large negative linear compressibility (NLC) was observed in metal-organic framework material Zn(HO3PC4H8PO3H)∙2H2O (ZAG-4) in experiment. We find a unique NLC mechanism

in ZAG-4 based on first-principle calculations. The key component to realize its large NLC is the deformation of H3O+ tetrahedron. With pressure increase, the oxygen apex approaches and then

is inserted into the tetrahedron base (hydrogen triangle). The tetrahedron base subsequently expands, which results in the b axis expansion. After that, the oxygen apex penetrates the

tetrahedron base and the b axis contracts. The negative and positive linear compressibility is well reproduced by the hexagonal model and ZAG-4 is the first MOFs evolving from non re-entrant

to re-entrant hexagon framework with pressure increase. This gives a new approach to explore and design NLC materials.

Most materials contract in all directions when hydrostatically compressed. That is the volume (), area () and linear () compressibility are all positive1. Negative volume compressibility is

thermodynamically impossible2. Counterintuitively, negative linear compressibility (NLC) indeed occurs in rare and remarkable crystals1,3,4,5,6. After screening of reported elastic constant

tensors from approximately five hundred crystals, Baughman et al. uncovered only 13 compounds showing negative compressibility in simple inorganic or organic compounds1. Among the thirteen

crystals, 11 structures were of monoclinic or lower symmetry. The typical positive linear compressibility (PLC) for crystal material lies in the range Kl ≈ 5–20 TPa−1, with lattice parameter

contracting 0.5~2% for each GPa increase in pressure7. Unfortunately, experimentally determined NLC, for a long time, had been below −2 TPa−1 (−0.2 TPa−1 for LaNbO48, −1.2 TPa−1 for Se9 and

−2 TPa−1 for BAsO410). Until recently, stronger NLC behavior is found: −3.8 TPa−1 for methanol monohydrate from 0 to 0.5 GPa at 160K11, −6.4 TPa−1 for α-BiB3O6 from 0 to 6.5 GPa12, −12 

TPa−1 for KMn[Ag(CN)2]3 from 0 to 2.2 GPa13, −41 TPa−1 for [Fe(dpp)2(NCS)2]∙py from 0 to 0.5 GPa14, −42 TPa−1 for Zn[Au(CN)2]2 from 0 to 1.8 GPa15 and −75 TPa−1 for Ag3[Co(CN)6] from 0 to

0.19 GPa7,16. Contrary to conventional materials, a specific direction of NLC material could not only increase with the increase of hydrostatic pressure, but also remain invariant15.

Therefore, NLC is a highly desirable property exploitable in the development of artificial muscles17, extremely sensitive pressure detectors, shock resistance materials and etc.

Metal-organic frameworks (MOFs) with extreme surface area and tunable pore structure have revolutionized the field of crystal engineering18,19,20,21. They consist of metal ions and organic

linkers, exhibiting various unique physical and chemical properties for diverse applications22,23,24. The crystalline order between metals and ligands combined with cooperative structural

transformability, forming flexible and responsive MOFs, namely soft porous frameworks. These materials can respond to mechanical stimuli in a tunable and precise fashion by molecular design,

which does not exist for other known solid-state materials25,26. The elastic behaviour of soft porous crystals is usually complex, such as anisotropic flexibility, negative Poisson’s ration

and high NLC25,26,27. The investigation of MOFs structure deformation under pressure can not only reveal the mechanism of these behaviours, but also can help us design new MOFs with desired

mechanical property. The research of the negative linear/area compressibility of framework materials started very recently28,29,30,31. The first case of NLC in MOFs was [NH4][Zn(HCOO)3],

which showed a high degree of mechanical anisotropy and negative compressibility Kl = −1.8 TPa−1 along its c axis from 0 to 0.94 GPa32. After that, NLC was found in silver(I)

2-methylimidazolate with Kl = −4.3 TPa−1 (along c axis, from 0 to 1 GPa)33 and [Ag(en)]NO3-I with Kl = −28.4 TPa−1 (along a axis, from 0 to 0.92 GPa)34. Clearfield and others pioneered named

the MOFs formed from the linker molecules with alkyl chains as zinc alkyl gate(ZAG) because of the likeness of the structure to a child safety gate (Fig. 1(a))21,35. Recently, Gagnon et al.

measured the lattice parameters of ZAG-4 under pressure with single crystal X-ray diffraction and observed NLC36,37. The b axis of ZAG-4 increases almost 2% in the range of 1.65–2.81 GPa,

indicating a strong NLC (Kl ≈ −16 TPa−1). Due to the inherently small atomic scattering factor of hydrogen, the exact positions of H2O in ZAG-4 can not be easily detected by X-ray

diffraction technique38. Aurelie U Ortiz et al. calculated the ZAG-4 structure under pressure, found a proton transfer and attributed the NLC (from 1.65 to 2.81 GPa) to this structural

transition39. However, PLC was observed after NLC in experiment. This explanation did not answer why NLC and PLC occurs subsequently after proton transfer. Therefore, an unambiguous

mechanism for the NLC in ZAG-4 is still an unresolved matter.

(a) ZAG-4 viewed along [0 0 1] direction. (b) ZAG-4 viewed along [0 1 0] direction. (c) Partial enlargement of (a), viewed approximately along [1 0 1] direction. At zero pressure, H1 is

close to PCO3 and far away from H2O. (d) H1 moves away from PCO3 with the increase of pressure and then H3O+ tetrahedron is formed. (e) Side view enlargement of H3O+ tetrahedron. (f) With

further pressure increase, apex oxygen penetrates H3O+ tetrahedron base and moves to the other side of hydrogen plane.

Density functional calculation, an integral part of MOFs research, is complementary to experimental techniques and offers invaluable information in characterization and understanding of

systems27,40,41,42. In order to elucidate the NLC mechanism of ZAG-4, we performed the density functional calculation using both PBE and Wu-Cohen (WC) functional43 as implemented in the

Quantum Espresso package44 to determine their atomic structures under pressure. The PBE functional, usually overestimating lattice parameters, has been widely used in density functional

calculations and the WC functional is known to be accurate in predicting solid volumes45,46,47. The wave function was expanded in a plane-wave basis set with an energy cutoff of 70 Ry and

the first Brillouin zone was sampled on a 3 × 3 × 4 mesh. The ultra-soft pseudopotential was used to represent the electron-ion interaction.

The conventional unit cell of ZAG-4 is depicted in Fig. 121,35,36. It has a base-centered monoclinic lattice with the b axis perpendicular to the a-c plane. Herein, the three building blocks

of ZAG-4 are the inorganic Zn-O-P-O chains and two bridging ligands (C4H8 and H2O). The 1D Zn-O-P-O chains orient along the c axis and are linked with each other along the b axis by H2O

molecules. So as to give a detailed description of the linkage between H2O molecules and Zn-O-P-O chains, we enlarge the bridging zone in Fig. 1(c). The ZnO4 and PO3C tetrahedrons are linked

along the c axis by sharing oxygen atoms and form the Zn-O-P-O chains. The water molecules are located between two Zn-O-P-O chains. As a result, Zn-O-P-O chains bridged by H2O form an

inorganic 2D structure parallel to the b-c plane. Furthermore, the inorganic planes are linked with each other along the a axis by C4H8, which is directly bonded with the P cations. In

brief, the 3D framework is established with the inorganic Zn-O-P-O chains extending along the c axis linked with each other by C4H8 chains and H2O molecules along the a and b axis,

respectively.

Firstly, we use the experimental crystal structure at zero pressure as our initial point and fully relax the lattice parameters and atomic positions with PBE functional. The calculated

lattice parameters (a = 19.00 Å, b = 8.41 Å, c = 8.18 Å and volume = 1214.27 Å3) agree well with experimental results (a = 18.51 Å, b = 8.29 Å, c = 8.27 Å and volume = 1160.55 Å3). As the

calculated volume is slight larger than the experimental value, our calculated bulk modulus of (11.6 GPa) is slightly lower than the experimental result (11.7 GPa). Although it is about one

thirty-fifth of the bulk modulus of sp3 carbon allotrope (around 400 GPa)48,49, this value is higher than that of porous MOFs MIL-53 and NH2-MIL-5330,50, but is lower than that of dense

MOFs51,52,53,54.

In order to explore its NLC mechanism, we applied hydrostatic pressure to ZAG-4 and investigate its structure variation. The calculated lattice parameters, accompanied with available

experimental data36, are shown in Fig. 2. We increase pressure and take the previously optimized structure as the initial point of higher pressure condition. In this way, we increase

pressure to 8 GPa and optimize the structure step by step.

Red and blue shaded areas in (b) manifest the NLC zone in experiment (from 1.65 to 2.81 GPa with increase of 1.8% in b axis) and calculation (from 5 to 6.25 GPa with the increment of 1.4% in

b axis), respectively.

As can be seen in Fig. 2, the experimental lattice parameters are well reproduced by our calculation based on the PBE functional. It is well-known that PBE calculation typically

overestimates the lattice parameters by 1–2%. As far as the overestimation is concerned, our calculated results agree excellently with the experimental data. As pointed by Gagnon et al., the

alkyl chains serve like a spring cushion and hence contract much under pressure36. With pressure increasing from zero to 2 GPa, the calculated volume compressibility is around 71 TPa−1 and

that of the experimental value from zero to 1.65 GPa is around 69 TPa−1. The calculated b and c axis have a jump from 2.25 to 2.5 GPa, which is accompanied with the proton transfer. Below

2.25 GPa, the H1 is close to the PCO3 octahedron and far away from the H2O (Fig. 1(c)). When pressure increases to 2.5 GPa, the proton turns to be close to the H2O and forms the H3O+

tetrahedron. Aurélie U Ortiz et al. also found this proton transfer and the H3O+ tetrahedron formation36,39. They used the wine-rack motif to explain the NLC of ZAG-4 after the proton

transition. However, as show in Fig. 2b, NLC does not occur immediately after proton transfer. Instead, the b axis smoothly expands from 5 to 6.25 GPa (blue shaded area in Fig. 2(b)), which

is far away from the proton transfer pressure (2.5 GPa). This indicates that proton transfer is not enough to lead to the NLC of ZAG-4. There must be something new.

We repeat these calculations with the WC functional, because the lattice constants of solids as determined by it are between LDA and PBE results and on average closer to experiment45,46,47.

In the WC results of ZAG-4, the H3O+ tetrahedron is formed at zero pressure. Consequently, the jump from 2.25 to 2.5 GPa in Fig. 2 does not exist in WC functional results (Supplementary

Information). Although the H3O+ is formed at zero pressure, the NLC does not take place immediately from zero pressure. Instead, the NLC of the b axis occurs in WC functional results from

3.5 to 4.75 GPa. This means the H3O+ tetrahedron and the NLC of b axis are reproduced in WC functional results as well.

We now pay our attention to the changes under pressure of the b axis. The experimentally observed NLC is from 1.65 to 2.81 GPa with the b axis increasing 1.8% (red shaded area in Fig.

2(b))36. The experimental data are too few (only two) to give a detailed description of the lattice parameters evolution. So as to draw a complete picture, we calculate the lattice

parameters from 5 to 7 GPa with a small pressure step of 0.25 GPa. As shown in Fig. 3(a), the b axis increases from 5 GPa and reaches its maximum at 6.25 GPa. Accordingly, the average NLC

from 5 to 6.25 GPa is −11 TPa−1, as strong as that of KMn[Ag(CN)2]3 (−12 TPa−1, from 0 to 2.2 GPa)13. After that, the b axis decreases with pressure increase and hence the linear

compressibility turns to be positive.

(a) Evolutions of the distance (dO-H triangle) between apex oxygen and H triangle and the b axis. (b) Area of the H triangle. With pressure increase, the apex oxygen approaches the H

triangle and expands its area and the b axis. After the apex oxygen penetrates the H triangle, the area of the H triangle decreases and the b axis contracts.

It is the deformation of the H3O+ tetrahedron that leads to the NLC of the b axis. As shown in Fig. 1(d,e), the apex oxygen of the H3O+ tetrahedron is above the H triangle at 5 GPa. The

distance (dO-H triangle) between the apex oxygen and the H triangle is defined to be positive at this condition (left inset of Fig. 3(a)). Comparing the left and right inset of Fig. 3(a) (or

Fig. 1(e,f)), we find that dO-H triangle turns from positive to negative with pressure increase. At the same time, evolution of dO-H triangle is completely accompanied with the evolution of

the b axis from expansion to contraction. We also calculate the structure of dehydrated ZAG-4. By deleting the H2O molecules in ZAG-4 at zero pressure, we get the initial structure of

dehydrated ZAG-4. Following the process of Fig. 2, we increase pressure and optimize the structure step by step. The calculated results (Supplementary Information) show that the b axis of

dehydrated ZAG-4 decreases smoothly with pressure increase. This means that the H3O+ tetrahedron deformation is essential to the NLC of ZAG-4.

It is easy to understand the NLC mechanism of ZAG-4 with the following picture in mind. Initially, dO-H triangle of the H3O+ tetrahedron is positive at low pressure. The apex oxygen moves

towards the H triangle with pressure increase. It leads to the decrease of dO-H triangle (Fig. 3(a)), the area expansion of hydrogen triangle (Fig. 3(b)) and the expansion of the b axis

(Fig. 3(a)). Therefore, the NLC of the b axis is directly results from the H3O+ tetrahedron deformation. With further pressure increase, the apex oxygen penetrates the H triangle and dO-H

triangle turns to be negative (right inset of Fig. 3(a)). From then on, the apex oxygen moves away from the H triangle. Consequently, both the area of the H triangle and the b axis decrease

with pressure increase. The calculated results using Wu-Cohen functional also show these characteristic features of H3O+ tetrahedron deformation. In brief, the b axis is expanded by the apex

oxygen approaching the H triangle and is contracted by the apex oxygen moving away from the H triangle.

There are four microscopic mechanisms frequently used: (i) ferroelastic phase transition, (ii) polyhedral tilt, (iii) helical system and (iv) wine-rack, honeycomb or related topology.6 The

NLC and PLC mechanism of ZAG-4 can be explained by the hexagonal model55,56. The non re-entrant and re-entrant hexagons are illustrated in the left and right inset of Fig. 4, respectively.

The linear compressibility Kl of this mode can be given analytically as56

The left inset is a non re-entrant hexagon with θ > 0, while the right one is a re-entrant hexagon with θ  θ > 3°) and the right side is PLC (θ