Synthesis of PCN–224

Micro-sized crystals of PCN–224 were prepared for TR-SFX analysis following a reported procedure38 with a slight modification. All of the chemicals required for the preparation of the MOF crystals were purchased from commercial suppliers (Sigma Aldrich and TCI Korea) and used as received without further purification. For the synthesis of the micro-sized crystals, ZrCl4 (0.21 g, 0.90 mmol; Sigma Aldrich, 99.5%), TCPP (0.070 g, 0.088 mmol; TCI Korea, 97.0%), benzoic acid (3.66 g, 30 mmol; Sigma Aldrich, 99.5%) and acetic acid (3.4 ml, 60 mmol; Sigma Aldrich, 99.7%) were suspended in N,N-dimethylformamide (DMF; 14 ml) in air in a 30 ml vial with a Teflon-lined cap and heated at 130 °C for 3 days to give red-brown cubic crystals. The vial was cooled for 1 h and the mother liquor was decanted. The remaining crystals were washed with DMF (10 × 10 ml) for 2 days and then with dichloromethane (DCM; 10 × 10 ml) for a further 2 days. To remove any solvents absorbed inside the pores of the resulting PCN–224, the washed crystals were heated at 150 °C for 12 h under vacuum. The characterization of PCN–224 is described in detail in Supplementary Methods.

Synthesis of PCN–224(FeII) and PCN–224(FeII)–CO

FeCl2 (558.0 mg, 4.4 mmol; Sigma Aldrich, 98%) and 2,6-lutidine (102.5 μl, 0.880 mmol; Sigma Aldrich, 97%) were mixed in DMF (5 ml) in a 25 ml round-bottomed flask. PCN–224 (150 mg, 0.037 mmol) was then immersed in this solution and heated for 12 h at 150 °C. After the reaction, the supernatant was decanted and the remaining crystals were washed in DMF (3 × 10 ml) for 30 min at 150 °C and then decanted. Next, 10 ml DMF was added to the washed crystals in the flask, which was then sealed and sonicated for 30 min. These crystals were washed with DMF (10 × 10 ml) and DCM (10 × 10 ml) for 2 days each. The remaining crystals were then dried for 12 h at 150 °C under vacuum, yielding dark-purple crystals. Unless otherwise noted, these procedures were conducted under N2 atmosphere using a Schlenk line or under Ar atmosphere in a glove box. After the activation of PCN–224(FeII), PCN–224(FeII) was exposed to 1 atm CO for 4 h to generate PCN–224(FeII)–CO.

Sample preparation for TR-SFX

Freshly activated PCN–224(FeII) was transported as a dry solid encapsulated in Ar to the Pohang Accelerator Laboratory X-ray free-electron laser (PAL-XFEL). A fixed-target sample holder was especially designed to keep the sample exposed to CO during the experiment. A fixed-target sample holder has two major advantages over conventional injector-style crystal delivery systems. First, the amount of sample consumed in the experiment can be greatly reduced. Second, there are no problems of chemical compatibility of the sample crystals with the medium because the crystals are simply placed on a substrate, such as a thin film, rather than dispersed in a liquid or gel-type medium. Acrylic Kapton tape was placed on one side of the sample holder and activated PCN–224(FeII) crystals were spread over the Kapton tape under an Ar atmosphere. Then, 50 μm polypropylene film was attached to the other side of the holder using double-sided acrylic tape. Samples encapsulated in the holder were exposed to ~1 atm CO for 4 h to generate PCN–224(FeII)–CO. After the exposure, the holder was completely sealed and subjected to TR-SFX analysis.

TR-SFX data collection

The diffraction experiments were performed at the XSS beamline of PAL-XFEL. The PAL-XFEL delivered X-ray pulses with a photon energy of 14.5 keV at a repetition rate of 30 Hz. The X-ray beam was focused to a spot of ~11 μm diameter full-width at half-maximum (FWHM). Laser pulses with a wavelength of 800 nm were generated from a Ti:sapphire regenerative amplifier and then converted to 400 nm pulses with a duration of <100 fs. The laser beam was focused to a spot of 45 μm diameter FWHM with a laser fluence of 5 mJ mm−2. Large-diameter X-ray and laser beams were used to prevent damage to the film and tapes and to keep the sample sealed in the CO environment. The X-ray and laser beams were aligned to overlap at the sample position with a crossing angle of 10°. The samples were delivered in home-made fixed-target sample holders mounted on Kohzu XA05A-L202 stages and moved for raster scanning. Time-resolved diffraction data were collected at 33 time delays between X-ray and laser pulses: −3.9, −2.4, −0.6, −0.2, 0.1, 0.4, 0.6, 0.85, 1.1, 1.35, 1.6, 1.85, 2.1, 2.35, 2.6, 2.85, 3.1, 3.35, 3.6, 4.1, 5.1, 7.1, 10.1, 17.9, 31.7, 56.3, 100, 178, 316, 562, 1,000 and 3,000 ps. More than 15,000 diffraction images were collected at each time delay. The diffraction images were recorded with a Rayonix MX225-HS detector in 4 × 4 binning mode (pixel size: 156 μm × 156 μm).

Data processing and analysis

We used the Cheetah software39 to select images containing meaningful diffraction spots by filtering out empty images and find spot positions in the selected diffraction images. Around 10–25% of the images collected from the sample at each time delay were identified as containing Bragg spots (Supplementary Tables 1–5).

The hit images obtained were further indexed and integrated using CrystFEL40. For indexing, the XGANDALF algorithm provided by CrystFEL was used, and the data beyond 0.9 Å were neglected due to the large experimental noise in the high-resolution region. The parameters for the experimental geometry were optimized by a trial-and-error approach. Intensities were integrated using the ‘rings-sat-cen’ option in CrystFEL with concentric rings of three, five and six pixels for the peak, buffer and background regions, respectively. The 'rings-sat-cen' integration option (i) utilizes three concentric rings for delineating peak and background areas, (ii) centers peak boxes on actual peak locations iteratively and (iii) does not integrate reflections that contain saturated pixels. The integrated intensities were merged using the partialator program from CrystFEL. XPREP was used to initially determine the space group of the crystal and prepare an input file for solving the structure. X-ray absorption was found to be negligible for crystals in the size range of 10–20 µm, where Tmin, the minimum value for experimental X-ray absorption correction, is larger than 0.9. Consequently, no absorption correction was applied. EDs and the corresponding molecular structures were solved using SHELXT41. VESTA42 was used to visualize ED maps overlain on refined molecular structures.

The space group of PCN–224(Fe)–CO is the same as that of a previously reported crystal structure of PCN–224(Fe) obtained by single-crystal X-ray crystallography38. The structure at each time delay, solved using SHELXT, was refined by a full-matrix least-squares method on F2 using SHELXL43. All non-hydrogen atoms were refined anisotropically and hydrogen atoms were refined isotropically. The structures were refined by applying displacement parameter restraints using DELU and SIMU commands and by applying a distance restraint using DFIX command in SHELXL. The disorder around the Zr clusters in the structures was modelled at each time delay. We observed a substantial residual ED in the ED maps that cannot be accounted for by the refined structures. This residual density may be attributed to either solvent molecules remaining in the pores or the partial occupation of MOF-525, which has a similar composition to PCN–224, but possesses a distinct structure. We employed the SQUEEZE option in PLATON44 to model the residual ED as contributions from disordered solvents or gas within the voids. The results of the refinement are summarized in Supplementary Tables 1–5. The A/B alerts generated from running CheckCIF and explanations are summarized in Supplementary Table 6. As the ORTEP-style drawings of the asymmetric unit of PCN–224(Fe)–CO at all time delays are very similar, we have included only the illustration corresponding to the laser-off condition as an example (Supplementary Fig. 1).

The experimental structure factor amplitudes (|F|) at each time delay were scaled based on the calculated (absolute) amplitudes. Δ|F| values were then obtained by subtracting the reference structure factor amplitude from the structure factor amplitudes at each time delay, represented as Δ|F|(hkl, t) = |F(hkl, t)| − |F(hkl, reference)|. Here, hkl represents the Miller indices and t denotes the time delay. The reference structure factor amplitude was derived from the average of the structure factor amplitudes measured at the negative time delays (−3.9, −2.4, −0.6 and −0.2 ps). DED (Δρ) maps were generated by Fourier transformation of each set of difference structure factor amplitudes using the phases from the PCN–224(Fe)–CO model before reaction45,46,47. Extended Data Figs. 1 and 2 show the DED maps at all time delays contoured at ±5.0σ (0.58 e Å−3) and ±3.0σ (0.35 e Å−3), respectively. Instead of a whole DED map, two regional DED maps, which are expected to arise from large structural changes, were used for further analysis: one was around the Fe atom and its CO ligand in the porphyrin ring (Extended Data Fig. 3) and the other was around the Zr6 node (Extended Data Fig. 4). In the DED maps obtained at late time delays, the original metal positions have strong negative densities surrounded by strong positive densities. Such features indicate that the metal atoms have been displaced from their original positions to other positions in random orientations, which is a characteristic of vibrationally hot structures. To confirm this hypothesis, the DED maps corresponding to the vibrationally hot structure, Ihot, were simulated for comparison with the experimentally observed DED maps. For the simulation, we calculated the ED maps for PCN–224(Fe)–CO before and after thermal excitation. For the ED map before thermal excitation, we used the ED map of the refined structure of laser-off PCN–224(Fe)–CO, which was obtained with the excitation laser off. To generate the ED maps of Zr and Fe atoms after thermal excitation, we increased the atomic displacement parameter of each atom in the refined structure of laser-off PCN–224(Fe)–CO and calculated the corresponding ED maps. The DED maps were obtained by subtracting the ED map before thermal excitation from that after thermal excitation (Extended Data Fig. 5). The simulated DED maps obtained on the assumption that the atomic displacement parameters increase due to the heating of the MOF have similar features to the experimental DED maps at late time delays (Extended Data Fig. 5). Based on these results, DED maps at late time delays were assigned to the heating component (that is, the vibrationally hot structure). Heat-free DED maps were obtained by removing the feature corresponding to the heating component from the DED map at each time delay (Supplementary Fig. 2).

Kinetic analysis of TR-SFX data

Data were collected at 33 time delays, which provides a large enough dataset for a quantitative kinetic and structural analysis. To obtain quantitative kinetic information from time-resolved DED maps, we applied SVD to decompose the original data into time-invariant DED features (that are, LSVs), their relative contributions (singular values) and their time profiles (that is, RSVs). First, the DED maps around the Zr6 node and the Fe porphyrin site were separately decomposed by SVD analysis to emphasize the dynamics of each metal domain. The singular values, autocorrelation values, and major LSVs and RSVs are shown in Extended Data Fig. 6. Inspection of the singular values and vectors revealed that, in both regions, the first component (the first LSVs and RSVs) stands out. The first RSVs for the two regions have almost identical shapes (Extended Data Fig. 7a). In addition to this dominant component, additional major components (the second RSVs) can also be observed in both regions.

For the primary component of both the Zr6 node and FeCO site, the laser-off positions of Zr and Fe exhibit strong negative densities surrounded by an isotropic shell of positive densities. This feature represents a vibrationally hot structure that extends through both metallic sites. Accordingly, the first RSVs represent the kinetics of the vibrationally hot structure, that is, thermal kinetics. To confirm this, we also performed a modified analysis using the Wilson plot method48. The Wilson plot was used to estimate the difference between the scale factors and the overall B-factors between datasets. The structure factors rescaled from 0 to 0.25 Å−2 during the structure refinement using SHELXL were used to extract the difference between the B-factors (Extended Data Fig. 7b). Specifically, we extracted the difference between the B-factor at each time delay and that of the averaged structure factors at negative time delays and compared the temporal trend with the first RSVs (Extended Data Fig. 7a). The first RSVs and the difference B-factors showed excellent agreement, indicating that the first RSV represents thermal kinetics. As the first RSVs for the two regions have nearly identical shapes (Extended Data Fig. 7a), they represent the kinetics of thermal excitation of the entire PCN–224(Fe) structure following photoexcitation. In this regard, we conducted a quantitative analysis of the time-resolved temperature changes in PCN–224(Fe); details of the analysis are provided in Supplementary Discussion and the results are presented in Supplementary Fig. 10.

For the Zr6 node, the second LSV exhibits an anisotropic feature along the d axis connecting the two different occupation sites of a Zr atom arising from positional disorder. Along the d axis, a negative density appears between the two occupation sites of a Zr atom, and two positive densities emerge at the exterior of the two occupation sites. This indicates that the Zr atoms were further disordered along the d axis. The components ranked third and below have low autocorrelation and singular values, and hence their contributions to the time-resolved DED maps around the Zr6 node can be ignored. For the FeCO site, the second LSV also shows clear anisotropic features along the axis connecting Fe and the CO ligand. Here, a pair of strong positive and negative densities appear above and below the original position of the Fe atom, indicating that this LSV represents the directional movement of the Fe atom. In addition, a strong negative density is evident around the position of the CO ligand, indicating that it dissociates upon photoexcitation.

The characteristic features of the RSVs can be briefly inspected before quantitative analysis. The second RSV for the Zr6 node shows three distinct behaviours: instantaneous rise, oscillation and exponential decay. As shown in Extended Data Fig. 6, the second RSV shows an instantaneous rise at around time zero, indicating that the contribution of the second LSV rises immediately upon photoexcitation. After the rise at time zero, the second RSV shows substantial oscillatory features with a peak at around 1.25 ps and a dip at around 4.0 ps. This oscillatory behaviour does not persist to late time delays, but damps quickly. Finally, the second RSV decays slowly towards zero. The second RSV of the FeCO site exhibits characteristics similar to those observed for the second RSV of the Zr6 node, particularly the instantaneous rise and exponential decay behaviours. However, the oscillatory behaviour of the second RSV of the FeCO site is notably weaker than that observed for the Zr6 node. Nevertheless, on closer examination, it is evident that these oscillatory features still make a discernible contribution. We performed a kinetic analysis for the two meaningful RSVs, the second RSVs of the Zr6 node and the FeCO site, to extract the detailed kinetics of the photoinduced anisotropic structural changes. For a quantitative analysis, the two second RSVs were globally fitted as a convolution of the IRF of ~200 fs FWHM with the sum of an exponential rise function (with a time constant of 20 fs, representing an instantaneous rise within the IRF of the experiment), an exponential decay function and a damped cosine function sharing time constants. The results of the fit are shown in Extended Data Fig. 8. A satisfactory fit was obtained with a time constant of 47.1 ± 0.5 ps for the exponential decay and two time constants of 5.55 ± 0.01 and 2.68 ± 0.02 ps for the period and damping, respectively, of the oscillatory motion.

Using the obtained time constants, we constructed a kinetic model that matches the behaviour of the major RSVs. The kinetic model consists of three species: (1) a species corresponding to the vibrationally hot structure (Ihot), (2) a species corresponding to the transient structure (Itr) that forms instantaneously within the IRF of the experiment upon photoexcitation and decays with a time constant of 47.1 ps, and (3) a species corresponding to the oscillatory motion (Iosc) with a period of 5.55 ps and a damping constant of 2.68 ps. The DED maps corresponding to these three species (that is, the SADED maps) were extracted via a kinetic analysis using the kinetic model. In this procedure, the entire DED maps of PCN–224(Fe)–CO were used instead of maps of just the regions around the Zr6 node and FeCO site. The kinetic analysis performed using the kinetic model is described in detail in Supplementary Methods. The resulting SADED maps for Iosc, Itr and Ihot and their time profiles are shown in Fig. 3a,b.

Generation of extrapolated maps

Detailed structural parameters for the structural species were obtained by structure refinement. The structure refinement was performed on diffraction intensities represented by the square of the structure factors of each structural species. The structure factors of the structural species, the so called extrapolated structure factors (Fextr(hkl))49, were obtained as follows: first, species-associated difference structure factors (Δ|F|SA(hkl)) for each structural species were obtained through inverse Fourier transformation of the corresponding SADEDs, derived from the kinetic analysis. The structure factors of the laser-off structure were used as the structure factors for ground state (FGround(hkl)). Then, the amplitudes of Fextr(hkl), |F|extr(hkl), were computed as follows:

$$|F_}}(hkl)=|F_}}(hkl)+(1/p)\times \,\varDelta |F_}}(hkl)$$

(1)

where |F|Ground(hkl) is the amplitude of FGround(hkl), p is the photoconversion yield ranging from 0 to 1 and Δ|F|SA(hkl) is the amplitude of ΔFSA(hkl). By determining the optimal value of p, |F|extr(hkl), free of the contribution of the ground state can be obtained. Suitable indicators for evaluating the level of photoconversion include the absence of negative density on the Zr atoms. As is evident from the term containing the division of ΔFSA(hkl) by p on the right-hand side of equation (1), it should be noted that when p is equal to unity (that is, when all of the ground state has been converted into the structural species), |F|extr(hkl) is simply the sum of the amplitudes of |F|Ground(hkl) and Δ|F|SA(hkl). When p is not equal to unity, the difference, Δ|F|SA(hkl), is weighted and added to |F|Ground(hkl) to obtain |F|extr(hkl). This process is referred to as extrapolation because typically, for p values less than 1, the small difference structure factors are extrapolated (or amplified to correspond to p = 1) and added to the structure factor of the ground state to determine the structure factors of the newly formed species. The p factors were determined to be 1.0 for Ihot, 0.20 for Itr and 0.25 for Iosc. Next, we calculated |F|extr(hkl) using the obtained |F|extr and the phase of the laser-off structure. The resulting |F|extr(hkl) values were squared to generate the diffraction intensities of the structural species. The diffraction intensities were refined in a similar manner to that used for the refinement of the laser-off structure, with additional restraints on the distances and planarity of the phenyl and pyrrole groups in the TCPP ligand. The A/B alerts associated with the structural species and the corresponding explanations are summarized in Supplementary Table 7.

The R factors for the extrapolated maps of the three structural species are comparable to, or even smaller than, those of the static crystal structures shown in Supplementary Tables 1–5. These small R factors for the structural species can be attributed to the use of a large number of DED maps for the calculation of the SADED maps. A SADED map can be described as a linear combination, or weighted average, of the DED maps, as detailed in Supplementary Methods. Due to this averaging process, in which the noise from each DED map is offset against one another, the noise in the SADED maps is reduced compared with that in individual DED maps. Consequently, SADED maps with reduced noise result in high-quality |F|extr, eventually leading to reduced R factors for the extrapolated maps. For instance, we can compare the quality of Δ|F|SA for Ihot and Δ|F| at 1 ns. At time delays much longer than the decay time constant of Itr (47.1 ps), only Ihot contributes to ΔF. Because p is equal to unity, it can be inferred that Δ|F| values at late time delays are approximately equivalent to Δ|F|SA for Ihot. In simple terms, the information from Δ|F| at multiple late time delays is averaged to extract Δ|F|SA for Ihot. Consequently, the quality of Δ|F|SA is inherently improved compared with Δ|F| for a single time delay, such as 1 ns. As Δ|F|SA of Ihot would be of higher quality than that of ΔF at 1 ns, it is reasonable to obtain a smaller R factor for the extrapolated map of Ihot (7.33%) than for the crystal structure at 1 ns (10.9%).