Nothing Special   »   [go: up one dir, main page]

$RTAWLOZ

Download as pdf or txt
Download as pdf or txt
You are on page 1of 40

International Journal of

Molecular Sciences

Review
Review of Applications of Density Functional Theory (DFT)
Quantum Mechanical Calculations to Study the High-Pressure
Polymorphs of Organic Crystalline Materials
Ewa Napiórkowska, Katarzyna Milcarz and Łukasz Szeleszczuk *

Department of Organic and Physical Chemistry, Faculty of Pharmacy, Medical University of Warsaw, Banacha 1,
02-093 Warsaw, Poland; enapiorkowska@wum.edu.pl (E.N.); kmilcarz@wum.edu.pl (K.M.)
* Correspondence: lukasz.szeleszczuk@wum.edu.pl

Abstract: Since its inception, chemistry has been predominated by the use of temperature to generate
or change materials, but applications of pressure of more than a few tens of atmospheres for such
purposes have been rarely observed. However, pressure is a very effective thermodynamic variable
that is increasingly used to generate new materials or alter the properties of existing ones. As
computational approaches designed to simulate the solid state are normally tuned using structural
data at ambient pressure, applying them to high-pressure issues is a highly challenging test of their
validity from a computational standpoint. However, the use of quantum chemical calculations,
typically at the level of density functional theory (DFT), has repeatedly been shown to be a great
tool that can be used to both predict properties that can be later confirmed by experimenters and to
explain, at the molecular level, the observations of high-pressure experiments. This article’s main goal
is to compile, analyze, and synthesize the findings of works addressing the use of DFT in the context
of molecular crystals subjected to high-pressure conditions in order to give a general overview of the
possibilities offered by these state-of-the-art calculations.

Keywords: pressure-induced phase transition; DFT calculations; polymorphism; solid state


Citation: Napiórkowska, E.; Milcarz,
K.; Szeleszczuk, Ł. Review of
Applications of Density Functional
Theory (DFT) Quantum Mechanical
Calculations to Study the
1. Introduction
High-Pressure Polymorphs of With decades of usage by physicists and mineralogists, high-pressure science is an
Organic Crystalline Materials. Int. J. established field. With the aid of improvements in experimental methods, high-pressure
Mol. Sci. 2023, 24, 14155. https:// research has expanded over the years from these roots into a variety of diverse domains.
doi.org/10.3390/ijms241814155 The field of crystal engineering is showing increasing interest in and research on the
Academic Editors: Qingwei Liao and behavior of organic and metal–organic systems under high pressure.
Xiuyu Wang The use of temperature to create or alter materials has dominated chemistry from
its very beginnings, while pressure of more than a few tens of atmospheres is only sel-
Received: 20 August 2023 dom encountered. However, pressure is a highly potent thermodynamic variable that is
Revised: 5 September 2023
increasingly applied to molecular systems to obtain new materials or modify existing ones.
Accepted: 11 September 2023
Depending on the topic of study, the term “high pressure” might refer to a number of
Published: 15 September 2023
different variables. On the scale employed by many planetary scientists, what is extraordinary
for biologists interested in life at large marine depths hardly qualifies. In the field of crystal
engineering, the range between approximately 0.1 and 20 GPa is usually described as “high
Copyright: © 2023 by the authors.
pressure”. It should be noted that the difference between atmospheric pressure and 1 GPa
Licensee MDPI, Basel, Switzerland. represents a factor of 104, which is an enormous value for molecular materials, where we can
This article is an open access article typically alter the temperature only within a few hundred degrees. Large structural changes in
distributed under the terms and molecular solids may be brought about by such pressure as contributions made by pressure-
conditions of the Creative Commons volume terms (pV) to free energy in this range can be equal to covalent bond energies [1].
Attribution (CC BY) license (https:// Numerous publications discuss the polymorphic changes caused by pressure in molecular
creativecommons.org/licenses/by/ solids. Sometimes, the high-pressure polymorph and the low-temperature form are identical,
4.0/). but this is not always the case. Moreover, high-pressure polymorphs may have the same space

Int. J. Mol. Sci. 2023, 24, 14155. https://doi.org/10.3390/ijms241814155 https://www.mdpi.com/journal/ijms


Int. J. Mol. Sci. 2023, 24, 14155 2 of 40

group symmetry as the original ambient-pressure form and even be isostructural with it, which
is defined as an isosymmetric phase transition. The only way to distinguish between a phase
transition and an anisotropic structural distortion is by direct optical inspection of the interface
and a series of structure refinements at various pressure values that demonstrate whether or not
the cell parameters and volume changes with pressure are continuous [2].
Since intermolecular potentials, atomic pseudopotentials, and other computational
protocols used to model the solid state are typically optimized with reference to ambient-
pressure structural data, their application to high-pressure problems is a very difficult
test of their validity from a computational point of view. Nevertheless, the application of
quantum chemical calculations, usually at the level of density functional theory (DFT), has
been proven multiple times to be an excellent tool that can be used not only to explain, at
the molecular level, the observations of high-pressure experiments [3], but also to predict
properties that can be further validated by experimentalists [4].
The main aim of this article is to gather, discuss, and summarize the results of the studies
addressing the application of DFT in the field of molecular crystals subjected to high-pressure
conditions to provide an overview of the possibilities offered by state-of-the-art calculations. Due
to the popularity and availability of DFT codes enabling such calculations, the complete list of
organic, inorganic, and metal–organic compounds that have been studied using computational
methods under high pressure is enormous. To describe all of these cases, a full-length book
would be required, which is beyond the scope of this review. Therefore, we limit the reviewed
works to those describing solid organics, following the title of this Special Issue.
This review is divided into sections. After the Introduction (Section 1), the main part
begins with a brief summary of the concepts of DFT calculations (Section 2), followed by
the presentation of the published results and a discussion of the most interesting, state-of-
the-art examples, grouped by the properties and applications of the studied compounds
(Section 3). Then, the following section presents the results dealing with fundamental
properties such as structure, enthalpy, vibrations, and phonons (Section 4), followed by a
part devoted to particular problems, e.g., anisotropic compression (Section 5), and finally a
section describing the general conclusions (Section 6).

2. Brief Theoretical Background


Atomistic model simulations based on DFT are undoubtedly the most popular quan-
tum chemical methods for the analysis of the structure–property relationships of crystalline
compounds [5–7]. Almost 60 years ago, by demonstrating that the total energy of electrons
moving in an external potential is a distinct functional of the electron density, Hohenberg and
Kohn (1964) [8] created the underlying formalism of DFT. Consequently, by decreasing the
total energy with regard to the electron density, the ground state of the quantum mechanical
system of electrons in a solid may be determined. Kohn and Sham (1965) [9] then proposed a
workable strategy for DFT computations, requiring the self-consistent solution of the so-called
Kohn–Sham (KS) equations. The KS equations translate an external potential system of n
interacting electrons into an effective potential system of n noninteracting electrons to describe
an external potential system of n interacting electrons. The Born–Oppenheimer approximation,
which isolates the dynamics of the electrons from those of the considerably heavier nuclei, is
used in the vast majority of implementations. DFT has been very popular for calculations in
solid-state physics since the 1970s. However, the rapid development and application of this
method in modeling molecular crystals started in the 1990s [10].
The type of functional and a basis set are the two most important variables that
must be chosen before starting the DFT calculations, and this choice greatly affects the
obtained results, as shown in Section 3.3.1 describing the studies of solid chlorothiazide.
Functionals characterize the electronic energy quantitatively, which, when combined
with the kinetic and electrostatic energy of the system, yields the system’s overall energy.
In studies devoted to organic solids, the most popular DFT functional is undoubtedly the
PBE [11] (Table 1). The choice of the basis set depends on the code being used. A range
of software basis sets (Table 2).
Int. J. Mol. Sci. 2023, 24, 14155 3 of 40

Table 1. Selected articles on the application of density functional theory (DFT) quantum mechanical calculations to study the high-pressure polymorphs of organic
crystalline materials.

Methods—DFT
Polymorphs Functional and
N◦ Molecule Pressure Range Type of Calculation Software Applied Year Ref. in Article
Studied Dispersion
Correction
A. Organic materials with metal additives
(IPy)4 (In2 Cl10 ) GO, HS, NBO, DOS, Crystal Explorer;
1. NP 0–1.51 GPa B3LYP 2021 [12]
IPy = 4-iodopyridinium OP, MO Gaussian
GO (compression
2. [Zn(µ-Cl)2(3,5-dichloropyridine)2 ]n P4b2, P4 0–9.34 GPa and decompression), CASTEP PBE TS 2021 [13]
Raman
Methylammonium lead bromide PBE-D3
3. Pm3m, Im3 0–2.5 GPa GO, BG Quantum ESPRESSO 2020 [14]
(MAPbBr3 ) Grimme
Tetragonal,
4. CH3 NH3 PbI3 (MAPbI3 ) orthorhombic, and 0–2 GPa BG, DOS VASP PBE 2019 [15]
cubic structures
Methylammonium lead bromide
5. Pm3m, R3m, R3 0–130 GPa GO, BG VASP PBE 2019 [16]
(MAPbBr3 )
Methylammonium lead bromide I (Pm3m) II (Im3), III
6. 0–3 GPa GO, BG, aiMD VASP PBE 2017 [17]
(MAPbBr3 ) (Im3), IV (Pnma)
Methylammonium lead iodide PBE-D3,
7. I4/mcm, Immm, Im3 0–1.95 GPa GO, BG NP 2016 [18]
(MAPbI3 ) B3PW91 + SOCPBE
8. Pt(bpy)Cl2, bpy = 2,20 -bipyridine Yellow and red form 0–3.8 GPa GO, MO, TD-DFT Gaussian B3LYP, LDA, BLYP 2007 [19]
B. High-energetic organic materials
USPEX, GO, BG,
9. 2,4,6-Trinitro-3-bromoanisole (TNBA) P21 /c, P21 21 21 0–10 GPa VASP PBE-D2 2022 [20]
DOS
Pentazolate anion (cyclo-N5 − ) salt
3,9-diamino-6,7-dihydro-5H-
10. bis([1,2,4]triazolo)[4,3-e:30,40- NP 0–50 GPa DOS, BG, PC, IR CASTEP PBE/G06 2022 [21]
g][1,2,4,5]tetrazepine-2,10-diium
((N5 − )2 DABTT2 + )
Pentazolate anion (cyclo-N5 − ) salt
11. N-carbamoylguanidinium NP 0–50 GPa DOS, BG, PC, IR CASTEP PBE/G06 2022 [21]
(N5 − GU+ )
Int. J. Mol. Sci. 2023, 24, 14155 4 of 40

Table 1. Cont.

Methods—DFT
Polymorphs Functional and
N◦ Molecule Pressure Range Type of Calculation Software Applied Year Ref. in Article
Studied Dispersion
Correction
1,3,5-Trinitrohexahydro-s-triazine
12. α, β, ε0 0–20.7 GPa GO VASP PBE vdW correction 2021 [22]
(RDX)
USPEX, XRD, GO, DIAMOND;
13. Ethylenediamine bisborane (EDAB) I, II, III 0–17 GPa vdW-DF, PBE 2021 [23]
PC, PF Quantum ESPRESSO
2,6-Diamino-3,5-dinitropyrazine-1-
NP; PBE;
14. oxide NP 0–25.7 GPa GO, DOS, BG AS, ES 2020 [24]
Gaussian B3LYP
(LLM-105)
15. 1,1-Diamino-2,2-dinitroethene (FOX-7) α, α0 , β, γ, δ, and ε 0–30 GPa GO, PC, Raman CASTEP PBE 2019 [25]
PBE-D2
16. 2,4,6-Trinitrotoluene (TNT) m-TNT and o-TNT 0–5 GPa GO CASTEP 2019 [26]
PBE TS
Octahydro-1,3,5,7-tetranitro-1,3,5,7- GO, HS AIM, IGM.
17. tetrazocine α, β 0–50 GPa MPD, MP, PC, PF, IR, CASTEP PBE 2019 [27]
(HMX) DOS
PBE-D2
18. Triaminotrinitrobenzene (TATB) NP 0–27 GPa GO, ZPE, Raman VASP 2017 [28]
Grimme
19. 1,1-Diamino-2,2-dinitroethene (FOX-7) α, α0 , ε 0–12.8 GPa GO CASTEP PBE Grimme 2016 [29]
20. Cyclotrimethylenetrinitramine (RDX) α, γ 0–10 GPa GO, TD, PC CP2K PBE-D3(BJ) 2016 [30]
2,6-Diamino-3,5-dinitropyrazine-1-
21. oxide NP 0–20 GPa GO, aiMD CP2K PBE-D2 2015 [31]
(LLM-105)
2,6-Diamino-3,5-dinitropyrazine-1-
GO, MO, aiMD (at 0 CASTEP;
22. oxide NP 0–45 GPa PBE-D2 2014 [32]
Gpa) CP2K
(LLM-105)
Octahydro-1,3,5,7-tetranitro-1,3,5,7-
B-HMX, insulator,
23. tetrazocine 0–130 GPa BG, QMD+ MSST CP2K SCC-DFTB 2014 [33]
metal
(HMX)
PWSCF; PBE,
GO, PC, TD, BG, CASTEP;
24. Silver fulminate (AgCNO) α, β 0–5 GPa PBE-D2; 2014 [34]
DOS WIEN2k
TB-mBJ
2,4,6-Trinitro-1,3,5-benzenetriamine PBE,
25. NP 0–7.02 GPa GO, CM Quantum ESPRESSO 2010 [35]
(TATB) PBE Grimme
Int. J. Mol. Sci. 2023, 24, 14155 5 of 40

Table 1. Cont.

Methods—DFT
Polymorphs Type of Functional and
N◦ Molecule Pressure Range Software Applied Year Ref. in Article
Studied Calculation Dispersion
Correction
Cyclotrimethylenetrinitramine 0–3.36 GPa and Quantum PBE,
26. α, γ GO, CM 2010 [35]
(RDX) 3.9–7.99 GPa ESPRESSO PBE Grimme
Hexanitrohexaazaisowurtzitane Quantum PBE,
27. β, γ, ε 0–2.7 GPa GO, CM 2010 [35]
(CL20, HNIW) ESPRESSO PBE Grimme
Quantum PBE,
28. Nitromethane (NM) NP 0–7.6 GPs GO, CM 2010 [35]
ESPRESSO PBE Grimme
Octahydro-1,3,5,7-tetranitro-
Quantum PBE,
29. 1,3,5,7-tetrazocine α, β, δ 0–7.47 GPa GO, CM 2010 [35]
ESPRESSO PBE Grimme
(HMX)
Quantum PBE,
30. Pentaerythritol tetranitrate (PETN) NP 0–9.04 GPa GO, CM 2010 [35]
ESPRESSO PBE Grimme
GO, PC, Raman, DMol3;
31. 1,3,5,7-Cyclooctatetraene (COT) NP 0 and 3.8 GPa PW91 2008 [36]
XRD, aiMD Gaussian; CPMD
Hexanitrohexaazaisowurtzitane α· H2O, β, γ, and CASTEP;
32. 0–400 GPa GO, SP, DOS, BG PBE; rPBE 2007 [37]
(CL-20, HNIW) ε DMol3
Supercell
approach
combined with the
33. Triclabendazole I, II 0–10 GPa embedded Gaussian ωB97XD 2022 [38]
fragment method,
GO, PC, IR,
Raman, TD
C. Pharmaceuticals
34. Chlorothiazide I, II 0–6.2 GPa GO, PC, TD, aiMD CASTEP PBE TS, PBESOL 2021 [39]
Quantum
35. Glycinium maleate NP 0–5.6 Gpa GO, PC, Raman LDA 2021 [40]
ESPRESSO
Int. J. Mol. Sci. 2023, 24, 14155 6 of 40

Table 1. Cont.

Methods—DFT
Polymorphs Type of Functional and
N◦ Molecule Pressure Range Software Applied Year Ref. in Article
Studied Calculation Dispersion
Correction
GO, DFT, and
DFTB3-D3(BJ) Quantum
36. Resorcinol α, β 0–4 GPa approach: ESPRESSO; B86bPBE-XDM 2021 [41]
vibrational Phonopy; DFTB+
frequencies
0–50 Quantum
37. Glycine α, β, γ GO, TD, BG PBE-D3 2020 [42]
GPa ESPRESSO
38. L-Histidine I, I0 , II, II0 0–7 GPa GO CASTEP PBE TS 2020 [43,44]
PBE TS, PBESOL,
39. Urea Form I and IV 0 and 3.1 GPa GO, PC, TD, aiMD CASTEP 2020 [45]
WC
USPEX, GO, ZPE, VASP, VASP VTST
40 A-Glycylglycine α, α0 , P21 21 21 0–18GPa PBE-D2 2020 [46]
PXRD tools for ZPE
41. L-Threonine I, I0 , II, III 0–22.31 GPa GO CASTEP PBE 2019 [47]
42. L-Threonine α, β 0–5 GPa GO CRYSTAL PBE-D3(BJ) 2019 [48]
PBE, PBE Grimme,
PBE TS, PBESOL,
PW91, PW91 OBS,
43. Glycine γ, δ 0–7.8 GPa GO, PC, TD CASTEP 2018 [49]
RPBE, WC,
CA-PZ, CA-PZ
OBS
GO (compression
VASP; PBE-D;
44. L-Serine I, II, III 0–8.2 GPa and 2017 [50]
Gaussian M06-2X
decompression)
Gaussian; M062X;
45. Tolazamide I, II 0–20 GPa GO, ZPE, TB 2017 [51]
VASP PBE-D3(BJ)
Quantum
GO, 2D PES, PC, B86Bpbe,
46. Aspirin I, II 0–12 GPa ESPRESSO, 2016 [52]
ZPE, TD B86bPBE-XDM
Phonopy
Int. J. Mol. Sci. 2023, 24, 14155 7 of 40

Table 1. Cont.

Methods—DFT
Polymorphs Type of Functional and
N◦ Molecule Pressure Range Software Applied Year Ref. in Article
Studied Calculation Dispersion
Correction
47. Aspirin I, II 0–5 GPa GO, PC, IR CRYSTAL B3LYP-2D 2015 [53]
48. Paracetamol I, II 0–5 GPa GO, PC, IR CRYSTAL B3LYP-2D 2015 [53]
WC, PBESOL,
GO, INS, PC, TD, CASTEP; PW91, PBE, rPBE,
49. Resorcinol α, β 0–4.5 GPa 2015 [54]
Raman CRYSTAL PBE-D2, PBE TS,
PBE/pob-TZVP
PBE, revPBE,
Quantum
50. Glycine α, β, γ, δ, ε 0–10 GPa GO vdW-DF, 2012 [55]
ESPRESSO
vdW-DF-c09x
PBE, revPBE,
Quantum
51. L-Alanine NP 0–10 GPa GO vdW-DF, 2012 [55]
ESPRESSO
vdW-DF-c09x
52. L-Serine I, II, III 0–8.1 GPa GO SIESTA PBE 2008 [56]
D. Others
53. Ammonium carbamate α, β 0–15 GPa GO CASTEP PBE TS 2022 [57]
54. Chloroform (CHCl3 ) P63 , Pnma 0–35 GPa GO, Raman, PF CASTEP PBE TS 2020 [58]
USPEX, GO
(compression and
55. Croconic acid Pca21 , Pbcm 0–55 GPa CASTEP PBE TS 2020 [59]
decompression),
PC, Raman, PF, BG
USPEX, GO
(compression and
56. Squaric acid P21 /m, I4m 0–25 GPa decompression), CASTEP PBE TS 2020 [59]
PC, Raman, PF,
BG, OP
Diisopropylammonium perchlorate PBE,
57. P1 0–3.3 GPa GO, Raman DMol3 2020 [60]
(DIPAP) PBE Grimme
Int. J. Mol. Sci. 2023, 24, 14155 8 of 40

Table 1. Cont.

Methods—DFT
Polymorphs Type of Functional and
N◦ Molecule Pressure Range Software Applied Year Ref. in Article
Studied Calculation Dispersion
Correction
g-C11 N4 , α-C11 N4 ,
58. C11 N4 d-C11 N4 , and 0 -70 GPa GO, PC, PF, TD VASP PBESOL 2019 [61]
β-C11 N4
Dihydrate, α and
59. Oxalic acid (−1.0)–12.0 GPa GO, XRD CASTEP PBE 2019 [62]
β
60. Sorbic acid C2/c 0–8 GPa GO, PC, Raman CASTEP PBE 2017 [63]
Oxobenzene-bridged Gaussian;
61. 1,2,3-bisdithiazolyl radical α, β, γ 0–13 GPa SP, MO, BG Quantum (U)B3LYP; PBE 2014 [64]
conductor (3a) ESPRESSO
[1a]2 : σ-dimer,
Bis-1,2,3-thiaselenazolyl radical
62. π-dimer 0-13.7 GPa GO, BG VASP PBE 2012 [65]
dimer
[1b]2 : NP
63. Indole HB and β 0–25 GPa GO, TB CASTEP PBE TS 2011 [66]
Abbreviations: 2D PES—two-dimensional potential energy surface; AIM—atoms in molecules; AS—absorption spectra calculations; aiMD—ab initio molecular dynamics; BG—bandgap;
CM—center-of-mass fractional position calculations; DOS—density of states; ES—excited state calculation; GO—geometry optimization; HS—Hirshfeld surface; IGM—intramolecular gradients
method; INS—inelastic neutron scattering; MO—molecular orbitals; MPD—mutual penetration distance; MSST—multi-scale shock technique; NA—not applicable; NBO—natural bond orbitals;
NP—not provided; OP—optical properties; PC—phonon DOS calculation; PD—phase diagram; PF—phonon frequency; QMD—quantum molecular dynamics; SCC-DFTB—self-consistent
charge density functional tight binding; SP—single-point calculations; SOC—spin–orbit coupling; TB—transition barrier calculation; TD—thermodynamics; TD-DFT—time-dependent density
functional theory calculations; VTST—variational transition-state theory; XRD—X-ray diffraction.
Int. J. Mol. Sci. 2023, 24, 14155 9 of 40

Table 2. Software applied in the reviewed papers in descending order of number of works.

N◦ Software/Code Basis Set Periodic License Type Ref. Method Number of Works Ref. in This Article
Academic, [13,21,25–27,29,32,34,37,39,43–
1. CASTEP Plane-wave 3D [67,68] 27
Commercial 45,47,49,54,57–59,62,63,66,69]
Academic,
2. VASP Plane-wave 3D [70,71] 11 [15–17,20,22,28,46,50,51,61,65]
Commercial
Free, General Public
3. Quantum ESPRESSO Plane-wave 3D [72,73] 9 [14,23,35,40–42,52,55,64]
License
4. Gaussian Gaussian-type orbitals Any Commercial [74,75] 8 [12,19,24,36,38,50,51,64]
Hybrid Gaussian-type Free, General Public
5. CP2K Any [76,77] 4 [30–33]
orbitals, plane-wave License
Academic,
6. CRYSTAL Gaussian-type orbitals Any [78,79] 3 [48,53,54]
Commercial
Numerically tabulated
7. DMol3 Any Commercial [80,81] 3 [36,37,60]
atom-centered orbitals
8. CPMD Plane-wave 3D Academic [82,83] 1 [36]
Slater-type orbitals,
Free, General Public
9. DFTB+ numerically tabulated Any [84,85] 1 [41]
License
atom-centered orbitals
Numerically tabulated Free, General Public
10. SIESTA 3D [86,87] 1 [56]
atom-centered orbitals License
FP-(L)APW + lo (the
full-potential (linearized)
11. WIEN2k 3D Commercial [88,89] 1 [34]
augmented plane-wave
and local orbitals)
Int. J. Mol. Sci. 2023, 24, 14155 10 of 40

Van der Waals interactions and hydrogen bonds, among other noncovalent forces,
are essential for the creation, stability, and operation of the majority of solid organics. At
present, high-level quantum chemical wavefunctions or the quantum Monte Carlo tech-
nique are the only ways to adequately account for omnipresent van der Waals interactions.
All prominent local-density or gradient-adjusted exchange-correlation functionals of DFT,
as well as the Hartree–Fock (HF) approximation, on the other hand, lack the appropriate
long-range interaction tail. The Kohn–Sham DFT equation does not take into account a
long-range electron correlation effect, sometimes referred to as the London component of
the dispersion energy factor. This problem has a long history of negatively impacting the
DFT calculations’ accuracy. There are several dispersion correction techniques accessible
today. They should be examined independently for each system of concern because their
inclusion does not always increase the calculation accuracy. Tkatchenko–Scheffler (TS) [90],
Grimme Dispersion (GD, often written as D) [91], and Many-Body Dispersion (MBD) [92]
are the three most used dispersion corrections in the study of solid organics.
The most frequent activity in all of the possible molecular simulations is the optimiza-
tion or minimization (with respect to potential energy) of the system being examined. The
ground-state structure of a compound is often calculated by minimizing the total energy
with regard to the locations of the atoms and, in the case of crystalline systems, also the
lattice parameters (a, b, c, α, β, γ) in athermal conditions, also omitting the zero point
motion (ZPVE) and at 0 GPa. This approach is usually called “full geometry optimization”.
The computation of structures at high pressure (but still within the athermal limit) is then
based on minimizing the enthalpy with respect to a non-zero stress tensor [93].
There are several methods to account for the impact of temperature in model cal-
culations, but each one needs a large amount of computing power. The free energy of
crystalline solids may be expressed using the quasi-harmonic approximation, in which
phonons do not interact with one another but the phonon frequencies depend on the unit
cell volume [94].
Molecular dynamics simulations are an alternate method to QHA that yield results that
go beyond the harmonic approximation. Conceptually, DFT-based molecular dynamics simu-
lations are identical to molecular-mechanics-based MD simulations in that the forces acting
on the atoms in a particular configuration are estimated, the atom locations are then updated
in accordance with the selected time step, and this process is continuously repeated [95].
The above paragraph only briefly summarizes the basic concepts of DFT calculations
and their practical implementations in describing molecular crystals, especially at high pres-
sure. For readers interested in more detailed information on the conceptual backgrounds of
such calculations, we recommend the more general, excellent review articles [3,93,96–99].

3. The Classes of Modeled Systems


The variety of modeled systems is enormous and there are multiple ways to group
them. In this review, we gather the materials into three main classes: (1) organic materials
with metal additives; (2) high-energy materials; and (3) pharmaceuticals. The first class
is on the border between inorganic materials, which are not the focus of this review, and
typical organic solids. Due to the presence of metals, the variety of bonds and interactions
existing in these materials is very large, which makes them interesting from a structural
point of view but also quite difficult to model. The second group is important in the context
of high pressure, which is usually connected with high-energy materials. The third group
has emerged as polymorphism that plays a key role in pharmaceutical sciences due to the
differences in stability, dissolution, and processability between the particular polymorphs
of active pharmaceutical ingredients.

3.1. Organic Materials with Metal Additives


Understanding the responses of organometallic materials to pressure changes is essen-
tial in optimizing their performance and designing new materials with tailored properties.
In silico methods provide a valuable tool for the investigation of the effects of pressure
Int. J. Mol. Sci. 2023, 24, 14155 11 of 40

on organic materials with metal additives, allowing for the efficient and cost-effective
exploration of a wide range of pressure conditions. The application of DFT methods has
appeared in many works as a useful tool in such cases (Table 1).

3.1.1. Methylammonium Lead Bromide (MAPbBr3 )


Among the group of organic materials with the addition of metals, MAPbBr3 has
attracted major attention. The DFT periodic calculations are known to be applicable in
describing the electronic properties of the structure. Therefore, the purpose of the first
research [14] was to investigate the structural flexibility of MAPbBr3 , which provides a
platform to engineer the optoelectronic properties by adjusting the bandgap energies by
applying external pressure or temperature. The band-edge calculations predicted the
bandgap to be direct up to 2 GPa, which is eligible in photovoltaics applications. According
to the results, the calculated bandgap energies were lower than the experimental values,
which is a well-known drawback of the PBE functional. Nevertheless, the absolute values
of the bandgap energies were not crucial for this analysis. More importantly, the reduction
of the bandgap with pressure observed in the DFT calculations was in agreement with the
experimental pressure-dependent photoluminescence energies in the cubic I (Pm3m) phase.
The authors of the other work [17] provided an experimental investigation of MAPbBr3
under various stress conditions supported by DFT calculations. First-principles molecular
dynamics (FPMD) simulations were performed to demonstrate the effect of the CH3 NH3
(MA) organic molecules on the instantaneous and time-averaged bandgaps of the various
structural phases. NVT FPMD simulations were carried out for pressure of 0.1/0.2/0.4 GPa,
0.6/1.0/1.4 GPa, and 1.8/2.4/3.0 GPa, on the phases Pm3m, Im3, Pnma, respectively. DFT
studies showed that the dynamics of the MA organic molecule and the inorganic lattice,
which were connected by N−H···Br hydrogen-bonding interactions, affected the distance
between the Pb and Br and the bandgap evolution under pressure. This study also high-
lighted the importance of MD simulations in obtaining good-quality results. According to
the static (at 0 K) calculation results, the bandgap of the Pnma phase at 0 K was predicted to
decrease as a result of pressure-induced band broadening. However, the bandgap obtained
by averaging the outcomes from snapshots of the FPMD simulations at 300 K increased
with increasing pressure, which was in agreement with the experimental observations.
The combination of experimental and theoretical investigations provided an insight into
the behavior of MAPbBr3 under pressure, which can be useful in the construction of
optoelectronic devices.
The aim of the next work [16] was to investigate band structure evolution on the
compression of halide perovskites like MAPbBr3 ; see Figure 1. According to calculations at
ambient pressure and confirmed experimentally, the investigated perovskites were direct
bandgap semiconductors. With the increase in pressure, bandgap contraction was observed
up to the first critical pressure at 4 GPa. Beyond this pressure value, the widening of
the bandgap happened and the bandgaps of organic perovskites became indirect. With
further compression, the bandgap remained direct; however, a decrease in the bandgap
was observed. The calculated results for MAPbBr3 were in good agreement with the
experimental findings, as was also the case in the earlier described works. With the help
of pressure, the magnitude and the direct–indirect nature of bandgaps can be adjusted. It
was pointed out that direct semiconductors are likely to have longer carrier lifetimes than
indirect ones, which is a desired property for materials in photovoltaic applications.

3.1.2. Methylammonium Lead Iodide (MAPbI3 )


In the next work [18], a similar material to the previously mentioned one was studied.
DFT calculations using hybrid functional B3PW91 with spin–orbit coupling (SOC) were
applied to compare the bandgap changes of MAPbI3 polymorphs under compression with
those obtained from laser-excited photoluminescence (PL) spectra. The calculated pressure-
dependent bandgaps decreased from 1.67 to 1.62 eV up to 0.3 GPa and then increased
from 1.58 to 1.74 eV between 0.4 and 2.0 GPa, which was in excellent agreement with the
Int. J. Mol. Sci. 2023, 24, 14155 12 of 40

experimental values. It was noticed that the transformation from I4/mcm to Im3 led to a
smaller bandgap energy, whereas conversion from the Im3 polymorph into the suggested
Immm polymorph reversed this tendency. Since photoluminescence is directly related to
polymorphism, conversion from the low-pressure I4/mcm phase to the cubic Im3 supercell
Int. J. Mol. Sci. 2023, 24, x FOR PEER REVIEW
(at 0.4 GPa) provided a distinct procedure for the reduction of MAPbI3 bandgaps and may
be useful in the development of materials in photovoltaic applications.

Figure 1. Band structures close to the Fermi level of cubic MAPbBr3 . Adapted with permission
Figure 1. Copyright
from [16]. Band structures close Chemical
2023 American to the Fermi level of cubic MAPbBr3. Adapted with permiss
Society.
[16]. Copyright 2023 American Chemical Society.
Another work [15] also focused on the structures of MAPbI3 under various external
pressure values. In the study, the structures of various phases of this material were op-
3.1.2.
timizedMethylammonium
and their electronic asLead
well Iodide (MAPbI
as optical 3)
properties were calculated. Calculations
included the density of states and partial density of states, band
In the next work [18], a similar material to the previously structuresmentioned
and projection
one was s
band structures, and the changes in charge transfer amount and electronic charge density.
DFT calculations using hybrid functional B3PW91 with spin–orbit coupling (SOC
At the phase transition point, the obtained bandgaps of these structures were higher than
applied to compare
their experimental the(error
values bandgaprange:changes of MAPbI
8%). Based polymorphs
on the 3results under
of this work, compressi
it was
those
possibleobtained from laser-excited
to evaluate changes photoluminescence
in optical properties (PL)
due to pressure-induced spectra.
phase The cal
transition.
pressure-dependent bandgaps decreased from 1.67 to 1.62 eV up to 0.3 GPa an
3.1.3. Chloroindium(III) Hybrid Perovskite (IPy)4 [In2 Cl]10
increased from 1.58 to 1.74 eV between 0.4 and 2.0 GPa, which was in excellent agr
In the next work [12], the authors investigated the structural and optical characterization
with the experimental values. It was noticed that the transformation from I4/mcm
of a chloroindium(III) hybrid perovskite (IPy)4 [In2 Cl]10 at high pressure. According to the
led to a smaller
experimental part ofbandgap energy,
the study, no whereaswere
phase transitions conversion
observed upfrom the
to 1.51 Im3
GPa. polymorph i
However,
suggested
several changes Immm polymorph
in electronic reversed
and structural this tendency.
properties were observed.Since
Basedphotoluminescence
on these findings, is
computational methods could help to understand the pressure effect
related to polymorphism, conversion from the low-pressure I4/mcm phase to the cu on the noncovalent
interactions and band structure, and thus the properties, of hybrid perovskites.
supercell (at 0.4 GPa) provided a distinct procedure for the reduction of MAPbI3 ba
and
3.1.4.may be useful
Zn(µ-Cl) in the development
2 (3,5-Dichloropyridine) 2 ]n of materials in photovoltaic applications.
Another
The authorswork of the[15]
nextalso
workfocused on thean
[13] performed structures
incremental ofsimulation
MAPbI3 under various e
of the hy-
pressure values. Inand
drostatic compression the study, the ofstructures
decompression of flexible
the plastically various phases of
coordination this materi
polymer
[Zn(µ-Cl) (3,5-dichloropyridine) ] . The tetragonal space group symmetry
optimized and their electronic as well as optical properties were calculated. Calcu
2 2 n was conserved
for anisotropic compression simulations, with the possibility of full cell relaxation, wherein
included the density of states and partial density of states, band structures and pro
an external potential was placed along the crystallographic c-axis. The anisotropic com-
band
pressionstructures,
was performed and the changes
to study how in thecharge transfer
compression amount
direction and the
affected electronic
observedcharge d
At thetransition.
phase phase transition
A series ofpoint, the obtained
hydrostatic compressionbandgaps of these
simulations structureswith
were performed were high
a pressure
their step of 0.5 GPa
experimental and no
values symmetry
(error range:restrictions.
8%). Based At around
on the 5 Gpa, a structural
results of this work,
possible to evaluate changes in optical properties due to pressure-induced
transition.

3.1.3. Chloroindium(III) Hybrid Perovskite (IPy)4[In2Cl]10


Int. J. Mol. Sci. 2023, 24, 14155 13 of 40

phase transition was observed, although the pressure-volume (p-V) curves did not show
discontinuities, typical for phase transformation. The DFT calculations revealed that this
transition was caused by the pressure-induced softening of low-frequency vibrations and it
was associated with breaking symmetry from P4b2 to P4. It was stated that the calculated
unit cell volume of the experimental ambient-pressure structure was overestimated by
less than 2%; thus, it was well represented by PBE TS. The results for compression and
decompression simulations were in good agreement with experimental data and they
showed that the phase transformation was predictable and reversible.

3.1.5. Pt(bpy)Cl2
The authors of the next article [19] applied DFT calculations to investigate the red-
to-yellow phase transition experimentally observed in Pt(bpy)Cl2 microcrystals at high
pressure. To simulate the yellow phase, the monomer was used, but to simulate the
stacked linear-chain red form, dimer and trimer units were employed. The theoretical
approach helped to analyze the experimentally observed differences in the Pt–Pt distance
and electronic structures manifested in differences in the optical spectroscopic properties of
the two phases. The presence of hysteresis in the optical properties allowed them to infer
that the yellow phase obtained on compression was not the same as observed at ambient
pressure. According to the calculations, the high-pressure luminescence corresponding to
the high-pressure yellow phase was more likely a halide–ligand transition rather than a
ligand–field transition. The authors concluded that supporting experimental methods with
theoretical work can be beneficial in the explanation of pressure-induced phase transitions.

3.2. High-Energetic Organic Materials


Studying high-energetic materials at various pressure values is crucial for several
reasons. Primarily, pressure has a significant impact on the stability, reactivity, and per-
formance of these materials. Detonation and deflagration expose energetic materials to
extreme conditions, resulting in vital transformations in their properties. These changes
impact various characteristics, including their susceptibility to shock initiation and their
chemical reactivity [100]. Additionally, studying high-energetic materials under different
pressure values can provide valuable insights into their structural changes, phase tran-
sitions, and energetic properties, enabling the development of safer and more efficient
energetic materials.

3.2.1. Silver Fulminate (AgCNO)


The authors of the next work [34] investigated the relative thermodynamic phase
stability of two polymorphs of silver fulminate (AgCNO) at high pressure and in a wide
range of temperatures (Figure 2). The research of energetic materials like AgCNO is
challenging due to their sensitivity to extreme conditions and risk of decomposition. Hence,
in silico simulations help to predict the physicochemical properties of these materials
under high pressure and temperatures. In this study, the authors used two plane wave
DFT codes, namely PWSCF and CASTEP, which is quite unusual as, in most of the other
works, the authors usually do not compare different programs. In addition, the authors
checked how the choice of the DFT functional and empirical dispersion correction affected
the results of calculations. This was beneficial as, according to the study’s results, the
pressure-induced phase transition from β to α at 2.5 GPa was observed only when the
geometry optimization was performed without any dispersion correction (GGA PBE).
These results were in contrast to those obtained using the DFT-D2 approach, indicating that
the α-form was more stable in the whole range of studied pressure values. Moreover, the
calculated volumes of the α and β forms were overestimated by 22.5% and 13.5% with GGA
PBE, for each phase, respectively. The application of DFT-D2 methods resulted in minor
discrepancies of approximately 2.2% and 4.5% for the α and β form, respectively. These
obtained results emphasize the need to include dispersion correction, which, in some cases,
is crucial to accurately model structural properties and their pressure dependence. Further,
These results were in contrast to those obtained using the DFT-D2 approach, indicating
that the α-form was more stable in the whole range of studied pressure values. Moreover,
the calculated volumes of the α and β forms were overestimated by 22.5% and 13.5% with
GGA PBE, for each phase, respectively. The application of DFT-D2 methods resulted in
Int. J. Mol. Sci. 2023, 24, 14155 minor discrepancies of approximately 2.2% and 4.5% for the α and β form, respectively. 14 of 40
These obtained results emphasize the need to include dispersion correction, which, in
some cases, is crucial to accurately model structural properties and their pressure
dependence. Further,
the bond length andthe bond
bond length
angle were and bond angle
evaluated were
up to evaluated
5 GPa up to 5 GPa
using DFT-D2 usingto
methods
DFT-D2 methods
characterize the to characterize the
compressibility compressibility
behaviors behaviors of the
of the polymorphs. polymorphs.
It was It was
noticed that the
noticed that the
polymorphs polymorphs
behaved behavedunder
anisotropically anisotropically
hydrostatic under hydrostatic
pressure. pressure.
The calculated The
electronic
calculated electronic
band structure usingband structure
the TB-mBJ using the
functional at TB-mBJ
ambientfunctional at ambient
pressure showed that pressure
the α and
showed
β phases were indirect bandgap insulators, with bandgap values of 3.51 and values
that the α and β phases were indirect bandgap insulators, with bandgap 4.43 eV,
ofrespectively.
3.51 and 4.43
MoreeV,examples
respectively.
of theMore examples
investigation of the investigation
of high-energy materialsofand
high-energy
the need to
materials andare
study them thedescribed
need to study
in thethem arebelow.
section described in the section below.

Calculatedenthalpy
Figure2.2.Calculated
Figure enthalpyasasa afunction
functionofofpressure
pressurefor
forα-α-and
andβ-polymorphic phasesofofAgCNO
β-polymorphicphases AgCNO
with (DFT-D2) and without (PBE-GGA) dispersion correction method. Reprinted from
with (DFT-D2) and without (PBE-GGA) dispersion correction method. Reprinted from [34], with the[34], with the
permission of AIP Publishing.
permission of AIP Publishing.

3.2.2. 3,5-Trinitrohexahydro-S-Triazine (RDX)


3.2.2. 3,5-Trinitrohexahydro-S-Triazine (RDX)
A few studies on RDX with the application of DFT methods have been published. One
A few studies on RDX with the application of DFT methods have been published.
of the works [69] investigated the compression behaviors of α, γ, and ε forms of RDX and
One of the works [69] investigated the compression behaviors of α, γ, and ε forms of RDX
the impact of pressure on the heat capacity of RDX using phonon calculations. However,
and the impact of pressure on the heat capacity of RDX using phonon calculations.
the results suggested very weak dependence in the case of all crystal forms of RDX,
However, the results suggested very weak dependence in the case of all crystal forms of
showing very small discontinuity at the α–γ phase transition. Furthermore, the authors
RDX, showing very small discontinuity at the α–γ phase transition. Furthermore, the
suggested that modeling pressure–heat capacity dependence is challenging and would
authors
requiresuggested that modeling
a very sensitive technique. pressure–heat
On the other capacity
hand, thedependence is challenging
obtained results and
of the lattice
would require a very sensitive technique. On the other hand, the obtained
parameters of α, γ were in good agreement with previous computational DFT-D studies results of the
lattice parameters of α, γ were in good agreement with previous computational
and experimental data. In this paper, the behavior of a high-temperature/high-pressure DFT-D
studies and experimental
polymorph, the ε-polymorph, data. was
In this
alsopaper, the behavior
analyzed of a high-temperature/high-
by a computational compression study
pressure polymorph,
for the first time and the theobtained
ε-polymorph, wasin also
results were good analyzed
accordanceby witha the
computational
experimental
compression
lattice parameters, with the largest deviation from the experimental volumeaccordance
study for the first time and the obtained results were in good being 0.9%.
with
Thethe experimental
geometry lattice parameters,
optimizations were performedwith the
by largest deviation
maintaining from the
the crystal experimental
space groups.
volume In being
another 0.9%.
workThe[22],geometry
the authors optimizations were performed
used DFT calculation by maintaining
to describe the differences thein
crystal space groups.
the intermolecular interactions and conformational changes of RDX under high pressure.
TheIn another work
combination [22], the authors
of experimental used DFT calculation
measurements to describe
and calculations the differences
indicated another phasein
the intermolecular
transition of β-RDX interactions
at 6.4 GPaand conformational
induced changes ofvariation
by the conformation RDX under high pressure.
of β-RDX.
The combination of experimental
In the next work measurements
[30], the authors calculatedand calculations
the Gibbs indicated
free energy another
to obtain phase
the relative
transition of β-RDX at 6.4 GPa induced by the conformation variation of
stability of forms α and γ of RDX at a pressure range of 0–10 GPa and temperature up to β-RDX.
450 In
K. the
Basednext
on work [30], the transition
the results, authors calculated the form
pressure from Gibbsαfree energy
to form γ was to suggested
obtain theto
relative stability of forms α andas,
be temperature-independent, γ ofinRDX at a pressure
the interval range
of range of 0–10
50–450 GPa
K, it wasand temperature
equal to 3.8 GPa,
which was in good accordance with the experimental outcomes.

3.2.3. 2,6-Diamino-3,5-Dinitropyrazine-1-Oxide (LLM-105)


Another commonly studied material is LLM-105. In the first work about this com-
pound, the authors used DFT calculations to support the experimental investigation of the
behavior of LLM-105, a nitro-energetic material, under high pressure [24]. Up to 25.7 GPa,
the first-order phase transformation was not observed. The DOS and band structure of
LLM-105 were calculated to determine the second-order phase transition at 10.5 GPa. The
combination of theoretical and experimental approaches yielded that an electronic struc-
Int. J. Mol. Sci. 2023, 24, 14155 15 of 40

ture phase transition occurred at ca. 10 GPa and was associated with a sudden change in
electronic transfer from the exposed oxygen atom to the amino groups. At around 25 GPa,
the abrupt bandgap reduction suggested a pressure-induced phase transition and this was
consistent with the experimental part of the work.
In the second work [31], it was also stated that up to 20 GPa, there was no indication
of a pressure-induced phase transition of LLM-105, even though the set of NPT molecular
dynamic simulations was performed at various pressure conditions to obtain the pressure
dependence on the lattice parameters. Nevertheless, at around 13 GPa, a significant shift
in the compressibility of the b-axis relative to the other axes was noticed, which was
attributable to a decrease in the spacing between molecules in adjacent in-plane sheets.
The authors concluded their work with a statement that EOS calculations can be used for
advancements in modeling high-pressure–temperature detonation reactions of energetic
materials such as LLM-105.
In the last study in this section [32], it was pointed out that the inclusion of dispersion
correction is crucial to accurately model systems like LLM-105. Geometry optimization of
the crystal structure at ambient pressure using the PW91 functional without dispersion
correction generated a ca. 27% error for the unit cell volume, mostly due to the larger errors
like 23.6% for the b-edge. However, using dispersion-corrected DFT methods allowed the
authors to obtain results in very good agreement with the experiments. The experimental
unit cell dimensions were well modeled up to 6 GPa, with a maximum deviation of
approximately 0.05 Å. Above 6 GPa, there were no available experimental data for this
moment to compare. There was no evidence of a structural phase transition up to 45 GPa
as well in this work. However, small symmetry changes and fluctuations in the structural
parameters were observed.

3.2.4. Cyclic Aliphatic Nitramine Octahydro-1,3,5,7-Tetranitro-1,3,5,7-Tetrazocin (HMX)


Another high-energetic material, HMX, has also obtained significant interest within the
literature [27]. The authors of one paper studied changes in the α-polymorph of HMX under
the influence of compression. They evaluated the pressure dependence on the hydrogen-
bonding network based on geometric evolution, electronic structure, Hirshfeld surfaces, the
method of atoms in molecules, the method of intramolecular gradients, the interpenetration
distance, mechanical properties, and vibrational spectra. The molecular structure, crystal
packing, lattice constants, volume, and bandgaps changed with increasing pressure, and
intramolecular interactions were dominated by hydrogen bonds. The results suggested
that the performance of the metastable α-HMX crystal is induced by pressure-triggered
hydrogen bond strengthening, and the intermolecular interactions clearly influence the
mechanical properties of the crystal. From phonon dispersion curves, a large negative
frequency at 4 GPa showed the dynamical instability of α-HMX at high pressure, which
was in line with experimental studies.
In another work [33], the pressure-induced metallization of condensed phase β-HMX
using quantum molecular dynamics in conjunction with a multi-scale shock technique
(MSST) was investigated and a self-consistent charge density functional tight binding (SCC-
DFTB) method was adopted. Based on simulations, a remarkable reduction in the bandgap
was observed with the increasing pressure exerted by the shock wave front, and the phase
became metallic at ca. 130 GPa, as a result of the breakage of the chemical bond, and
chemical chain reactions were initiated. The major gaseous products of the decomposition
of β-HMX were in agreement with experimental and previous theoretical results.

3.2.5. 1,1-Diamino-2,2-Dinitroethene (FOX-7)


Another high-energy material appearing in the literature is FOX-7 [25]. The authors
of one of the studies examined FOX-7 using dispersion-corrected density functional theory
under uniaxial compression along three major crystalline directions. Calculations included
principal and shear stresses as well as Raman spectra prediction to explain the role of hydrogen
bonding in stability and phase transitions under pressure. The simulated Raman spectra
Int. J. Mol. Sci. 2023, 24, 14155 16 of 40

showed anomalous changes attributed to the rearrangement of inter- and intramolecular


structures, the alternation of bond lengths, and the rotation of groups. The findings of this
work showed that non-hydrostatic pressure in comparison to hydrostatic compression could
trigger a phase transition, possibly along different crystallographic directions. This work
shows promising insights into studies performed on energetic materials such as FOX-7 in the
case of hydrogen bonding, anisotropy, and sensitivity to shock initiation.
In the second study [29], the authors investigated the structural response of an in-
sensitive energetic crystal of FOX-7 up to 12.8 GPa using synchrotron single-crystal X-ray
diffraction measurements and DFT calculations (Figure 3). The comparison of the ex-
perimental unit volumes and angles with theoretical findings showed that they were
reproduced well with the DFT-D approach, within a 1.3% error for the investigated pres-
sure range. The lattice parameters were also well modeled, both below and above 4.5 GPa,
with abrupt changes at the transition pressure. It was noticed that the predicted length
along the axis normal to the molecular layers (precisely the b-axis below 4.5 GPa and a-axis
above 4.5 GPa) varied from the measured values more than the lengths along the other two
axes. Based on these findings, it was confirmed that even using the DFT-D approach did
not reproduce well the weak interlayer interactions as well as strong in-layer interactions.
Due to the increase in strength interactions at higher pressure, the agreement between
the experimental values improved significantly. Overall, the results demonstrated that
structural changes in the FOX-7 structure were reproduced well by DFT-D calculations
of various phases and at the applied pressure range. The authors demonstrated that the
chemical and structural stability of FOX-7 was mainly controlled by the accommodation
of external compression through the anisotropic compressibility of the unit cell, phase
Int. J. Mol. Sci. 2023, 24, x FOR PEERtransformation
REVIEW of the structure to planar layers, and an increase in the hydrogen and 17 of 4
C−NO2 bonding strength. The authors stated that the results can be beneficial for a further
understanding of the shock responses of insensitive, highly explosive crystals.

Figure3. 3.High-pressure
Figure High-pressure structural
structural response
response of an of an insensitive
insensitive energetic
energetic crystal—1,1-diamino-2,2
crystal—1,1-diamino-2,2-
dinitroethene
dinitroethene (FOX-7);
(FOX-7); a, b,
a, b, c—unit
c—unit cellcell lengths.
lengths. Adapted
Adapted withwith permission
permission fromfrom
[29]. [29]. Copyright 202
Copyright
American
2023 AmericanChemical
ChemicalSociety.
Society.

3.2.6. 2,4,6-Trinitrotoluene (TNT)


Konar et al. studied two polymorphs of 2,4,6-trinitrotoluene (TNT), namely m-TN
and o-TNT, with the use of the PBE functional and two dispersion corrections, PBE TS an
PBE-D2 [26], to evaluate the agreement between the experimental and correspondin
Int. J. Mol. Sci. 2023, 24, 14155 17 of 40

3.2.6. 2,4,6-Trinitrotoluene (TNT)


Konar et al. studied two polymorphs of 2,4,6-trinitrotoluene (TNT), namely m-TNT
and o-TNT, with the use of the PBE functional and two dispersion corrections, PBE TS
and PBE-D2 [26], to evaluate the agreement between the experimental and corresponding
calculated results. The results showed that the PBE-D2 method outperformed PBE TS
in the case of a negligible temperature effect. At ambient conditions, though, PBE TS
showed an advantage over PBE-D2. These results suggested that the choice of dispersion
correction may depend on the specific pressure and temperature conditions of the system
being studied.

3.2.7. Pentazolates
In the next work [21], DFT calculations were used to investigate the geometry and crys-
tal structures, electronic features, hydrogen-bonding network, and vibrational properties of
two energetic pentazolate anion salts up to 50 GPa. Studying crystal structures’ behavior
and their properties under extreme conditions provides guidelines for the application of
pentazolates as energetic materials. Based on the results, in the case of (N5 − )2 DABTT2 + ,
the critical point was found at 9 GPa, indicating the possibility of a phase transition at this
pressure value. Furthermore, for N5 GU+ , it was not found. The observed distortions of the
cations of both compounds were investigated, as they could affect the hydrogen-bonding
network in the crystals and can have an impact on their stability. From the gradually de-
creasing bandgap under increasing pressure, it was noticed that the pressure enhanced the
charge overlap and delocalization, thus improving the electronic transition from occupied
orbitals to empty orbitals.

3.2.8. 2,4,6-Trinitro-3-Bromoanisole (TNBA)


In another work [20], the authors applied the first-principles evolutionary crystal
structure prediction method USPEX to evaluate the possible polymorphs of 2,4,6-trinitro-3-
bromoanisole (TNBA), an energetic material, at 10 GPa. In order to find the most stable
structures, the search was performed without any lattice or symmetry constraints, allowing
the structure to be flexible. Although the validation search was performed for a structure
stable at 0 GPa, this ambient phase crystal structure was not found. However, the structures
that were obtained at 0 GPa were energetically competitive, with the difference from the
ambient pressure phase being lower than 1 kcal/mol. An unreacted equation of state (EOS)
for the ambient phase of TNBA was found; with the optimization of the lattice parameters,
atomic configurations, and hydrostatic pressure, an isothermal method (at T = 0 K) was
performed and the result was presented in the form of the function of the volumetric com-
pression ratio V/V0 . Calculations were performed for equilibrium volume V0 at ambient
pressure and the reduction in the volume was sequential by 0.02 V0 . The authors stated
that the crystal structure prediction (CSP) simulations concerned decomposition products;
however, they did not calculate the energy barrier, which is important in understanding
the chemistry of energetic materials.

3.2.9. Hexanitrohexaazaisowurtzitane (HNIW or CL-20)


Another article focused on the polymorphs of CL-20 [37]. The β and ε-CL-20 unit cells
were optimized in CASTEP. However, due to the large size of α· H2 O phase crystals, it was
beyond the computational capabilities of the authors to perform its geometry optimization
at that time (in 2007) and the γ-phase could not be optimized for unknown reasons. Hence,
the single-point calculations of four forms were performed using the DMol3 software.
Nevertheless, based on bandgap calculations, the experimental order of phases was accu-
rately predicted. The least sensitive form (ε-CL-20) was used for further analysis under
hydrostatic pressure. At 0 GPa, the calculated lattice constants were 1.5% or so larger
than the measured ones, which is typical for the GGA-PBE method. At low and high
pressure, the compressibility of ε-CL-20 is anisotropic; up to 10 GPa, minor changes in the
lattice parameters, band structure, and DOS were observed. However, from 50 to 400 GPa,
Int. J. Mol. Sci. 2023, 24, 14155 18 of 40

they changed greatly. The increment in pressure reduced the bandgap and consequently
increased the sensitivity of ε-CL-20. Based on these changes, 400 GPa was considered to be
the critical pressure for the insulator–metal phase transition.

3.2.10. Triaminotrinitrobenzene (TATB)


In the next study [28], the authors applied DFT calculations to simulate Raman-active
modes and compared theoretical models with experimental spectra of TATB. The simulated
crystal structure at ambient pressure, as well as its volume–pressure dependence and Raman
spectra, using DFT methods with implemented dispersion correction and thermal and zero-
point energy contributions, were in good agreement with experimental findings. Furthermore,
the calculated peak positions were in excellent agreement with experimental measurements
from both this study and those performed by other researchers [101], with the highest inaccu-
racy of 10% in the low frequencies and less than 3% elsewhere. Up to 27 GPa, no first-order
transition was observed as no discontinuities within the low-frequency mode evolution were
detected, which was in accordance with the experimental study [28].

3.2.11. RDX, HMX, CL-20, NM, TATB, and PETN


Most of the earlier mentioned energetic materials were also studied in the next
work [35]. The authors performed a theoretical investigation using DFT calculations
at ambient pressure and on compression and compared the calculated results with the
experimental findings. In addition, sets of calculations were performed to test the impact
of applying dispersion correction and a selected energy cut-off from 25 to 80 Ry on the
accuracy of the predicted lattice parameters at ambient pressure. It was also emphasized
that the level of agreement was highly improved by using the DFT-D method instead of
the conventional DFT. Overall, good agreement between the calculated and experimental
structures’ parameters was obtained using dispersion-corrected calculations as the results
were within a 2% error range. Afterward, the hydrostatic compression effects on α and
γ-RDX, β-HMX, ε-CL20, NM, TATB, and PETN crystals were explored at a wide pressure
range selected to match the available experimental data, e.g., a range of 0–3.36 GPa for
α-RDX and 3.9–7.99 GPa for γ-RDX. Similarly, in this case, including dispersion correction
in the calculations improved the obtained results. On this basis, the authors claimed that
the DFT-D method was capable of reproducing the changes in lattice parameters and unit
cell volume caused by hydrostatic compression as the obtained maximum errors for the
lattice parameters were 1.8% (α -RDX), 1.07% (γ-RDX), 3.67% (β-HMX), 0.91% (ε-HNIW),
2.62% (NM), 1.74% (TATB), and 2.79% (PETN), respectively, relative to compression data
obtained at ambient temperature.

3.3. Pharmaceuticals
3.3.1. Chlorothiazide
It has been noticed that phase transition can be kinetically hindered and applying
pressure in static calculations may not be sufficient to detect it, even if the relative stability of
polymorphs could suggest that the phase transformation should occur at the investigated
conditions. To overcome this issue, several methods have been applied in theoretical
approaches. It is worth mentioning that these approaches are usually very computationally
demanding and their application can be limited by the size of the crystal’s unit cells to
small or moderate-sized systems.
One of the solutions to the aforementioned problem can be the application of ab-initio
molecular dynamics (aiMD) to include the kinetic energy factor in the phase transition [39].
For example, based on experimental findings, chlorothiazide, a diuretic agent, has been
proven to undergo pressure-induced isosymmetric structural phase transition (IPT) of
Form I to Form II at 4.2 GPa (Figure 4) [102]. The authors of the work [39] performed
geometry optimization calculations of Form I and Form II at all experimentally studied
pressure conditions to simulate compression (when geometry optimization starts from
Form I) or decompression (when geometry optimization starts from Form II) using two
Int. J. Mol. Sci. 2023, 24, 14155 19 of 40

functionals: PBE TS and PBESOL. The next stage was to measure the root mean square
deviation (RMSD) values between the estimated and experimental crystal structures in
order to assess the correctness of the modeled structures. According to the energy results
obtained using the PBE TS functional, Form II ought to be more stable at higher pressure.
However, the energy barrier prevented the IPT from being observed during geometry
optimization. In the case of the PBESOL, Form II was preserved only for pressure higher
than 3.5 GPa, and a jump discontinuity in the lattice parameters would suggest that, at this
point, one can expect IPT. However, according to the PBESOL calculations, Form I should
be more energetically favorable than Form II in the whole studied pressure range, which
is in contrast to experimental observations and indicates that the dispersion correction is
crucial for the accurate modeling of solid-phase transition. Hence, only the PBE TS was
used for the thermodynamic and aiMD calculations. According to the variations in the
predicted thermodynamic parameters (∆H and TS), the pressure-induced IPT should be an
entropy-driven transition. According to the computed changes in the Gibbs free energy
(∆G), the investigated forms should coexist at pressures between 3.5 and 4.1 GPa and,
Int. J. Mol. Sci. 2023, 24, x FOR PEER REVIEW 20 of 41
above this pressure, Form II should be more stable and dominant, which is in agreement
with the experimental data. Therefore, by adding kinetic energy in the aiMD simulations, it
was possible for the studied system (Form I) to overcome this energy barrier and reach the
static
deepercalculations,
minimum (Formthe authors successfully
II). This approach isapplied aiMD to overcome
very computationally energybut
demanding barriers
allowsand
observe the IPT of Form I to
us to overcome kinetic barriers. Form I during the simulation.

Figure
Figure 4.
4. Geometry optimization
Geometry optimization of two
of two polymorphic
polymorphic forms
forms (I (Iofand
and II) II) of chlorothiazide
chlorothiazide at various
at various pressure
pressure
using two DFT functionals (PBESOL and PBE + TS). Adapted from [39], licensed under CC BY 4.0.under CC
using two DFT functionals (PBESOL and PBE + TS). Adapted from [39], licensed
BY 4.0.
To conclude, the geometry optimization allowed the researchers to observe the pressure-
induced IPT of Form II to Form I during decompression. The thermodynamic parameters
3.3.2. Urea
indicated the possibility of reversed transformation and determined the relative stability of
In the next work [45], the application of aiMD was described in the case of urea
(Figure 5). The authors chose two approaches to predict the phase transition of urea
between Form I and IV, occurring at 3.1 GPa. Firstly, they performed calculations of
geometry optimization and thermodynamic calculations at 0 and 3.1 GPa. They stated that
the accuracy of the data depended on the chosen functional. At lower pressure,
Int. J. Mol. Sci. 2023, 24, 14155 20 of 40

the forms at higher pressure. Based on the promising results of the static calculations, the
authors successfully applied aiMD to overcome energy barriers and observe the IPT of Form I
to Form I during the simulation.

3.3.2. Urea
In the next work [45], the application of aiMD was described in the case of urea (Figure 5).
The authors chose two approaches to predict the phase transition of urea between Form I
and IV, occurring at 3.1 GPa. Firstly, they performed calculations of geometry optimization
and thermodynamic calculations at 0 and 3.1 GPa. They stated that the accuracy of the data
depended on the chosen functional. At lower pressure, calculations performed on Form IV
as the initial structure fit the results obtained when the initial form was Form I, with good
accuracy. A change in the crystal space group was observed as in the experiment, as well as
the final unit cell dimensions. The values of energy and free energy, which were obtained
in initializing calculations from different structures, were also almost identical. However,
geometry optimization performed at higher pressure did not result in a change in the space
group of Form I and the final structure differed significantly from the one that was obtained
when starting from Form IV, but the values of differences in the energy and Helmholtz
free energy between both structures correctly predicted their stability. The second approach
involved the application of NPT quantum molecular dynamic calculations at 3.1 GPa. The
received unit cell lengths of Form I equaled those for Form IV, as well as the final space group.
Int. J. Mol. Sci. 2023, 24, x FOR PEER REVIEW
The application of QMD calculations was stated to be a successful method in the study of the
phenomenon of polymorphic phase transitions and a useful tool for accurate crystal structure
prediction at different pressure values.

Figure 5. Ab initio molecular dynamics simulations of urea Form I at 3.10 GPa. Polymorphic phase
Figure 5. Ab initio molecular dynamics simulations of urea Form I at 3.10 GPa. Polymorph
transition is observed after 6 ps of simulation. Adapted from [45], licensed under CC BY 4.0.
transition is observed after 6 ps of simulation. Adapted from [45], licensed under CC BY 4.0
3.3.3. Tolazamide
3.3.3.The
Tolazamide
polymorphic transition of tolazamide, an anti-diabetic drug, is another example of
when a solid-state transformation is kinetically hindered. The authors of [51] investigated
The polymorphic transition of tolazamide, an anti-diabetic drug, is another e
the pressure dependence on two forms of tolazamide (Form I and II) and the transformation
of
to awhen a solid-state
denser form transformation
on compression is kinetically
using DFT calculations hindered.
and experimental The The
methods. authors
investigated the pressure
transition of polymorph I into dependence
Form II was noton two forms
observed of tolazamide
with increasing (Form
pressure, I and II)
neither
experimentally nor computationally. Furthermore, Form II did not transform
transformation to a denser form on compression using DFT calculations and exper into Form
I, even though Form I was determined to be denser at all pressure values discussed in
methods. The transition of polymorph I into Form II was not observed with inc
pressure, neither experimentally nor computationally. Furthermore, Form II d
transform into Form I, even though Form I was determined to be denser at all p
values discussed in this work. The authors proposed two hypotheses to explai
Int. J. Mol. Sci. 2023, 24, 14155 21 of 40

this work. The authors proposed two hypotheses to explain these observations. First,
there was no thermodynamic-driven force due to the similarity in the free energies of
polymorphs at the studied pressure. Second, Form I was more thermodynamically stable;
however, the solid-phase transition was kinetically controlled. To test the former, the
geometry optimization of two polymorphs was performed up to 20 Gpa. From the results
of calculations, Form I had lower internal crystal energy and enthalpy than Form II at all
simulated pressure values. The estimated enthalpy difference between the two phases was
6.1 kJ/mol and 96.8 kJ/mol at ambient pressure and 20 GPa, respectively. Interestingly,
the zero-point energy (ZPE) contribution to the difference between the enthalpies did not
exceed 4 kJ/mol and thus did not affect the obtained results. The discrepancy between the
calculated and experimental outcomes could have occurred due to the lack of an entropy
contribution in the calculated values. Nevertheless, the calculated results showed that Form
I was more stable across the entire studied pressure range; thus, no transformation of Form
I to Form II should be expected. On the other hand, the transformation of Form II to Form
I could be possible if there was no kinetic barrier. Further investigation of the pressure
impact on the structures of tolazamide revealed that the intermolecular interactions changed
significantly under pressure. Comparison of the conformations of the polymorphs showed
that the main difference was associated with the change in the C4–S1–C1–C8 dihedral angle;
thus, tolazamide was assumed to be a conformational polymorph. The calculated value
of the energy barrier of the transformation of one conformer into another was 28 kJ/mol,
and this was in agreement with previously reported data. In a solid state, the estimated
energy barrier could be even higher due to the intermolecular changes required to obtain
conformational changes. Under ambient conditions, the barrier can easily be overcome
in solutions, which was confirmed by the observation of the recrystallization of Form
II into Form I in dioxane and methanol. Additionally, it was noted that increasing the
temperature at ambient pressure facilitated the kinetically hindered phase transformation
of the tolazamide polymorphs.

3.3.4. Aspirin
The calculation of the energy barrier can help to justify the lack of experimentally
observed polymorphic transition under pressure, as in the case of aspirin crystals. The
behavior of these crystals under pressure was the area of interest in the next work [52]. The
authors used dispersion-corrected DFT calculations to establish the relative stabilities of two
forms of aspirin (Form I and Form II) and the shear–slip mechanism of their interconversion.
Based on the result, Form I was predicted to be 0.3 kJ/mol per molecule less stable than
Form II with B86bPBE-XDM. Including ZPE slightly decreased the difference by 0.1 kJ/mol,
which was in reasonable agreement with previous studies using the PBE TS and PBE MBD
approaches [103]. The predicted free-energy (∆G) difference of 0.3 kJ/mol per molecule at
298 K, as thermal effects were included, favored Form I as being more stable than Form
II, corresponding to the experimental preference. From the results of the PES scan for
the phase transition from Form II to Form I, the obtained energy barrier of the transition
increased from 10 kJ/mol per molecule at 0 GPa to the maximum value of ca. 18 kJ/mol
per molecule at 7 GPa [52]. The results explained the lack of transition from Form II to I,
supported by experimental work up to 10 GPa within limited timescales [104].

3.3.5. Triclabendazole
As a pressure increase can be insufficient to observe a phase transition, in some cases,
thermodynamics calculations should be performed to define the conditions required for
the occurrence of phase transformations [30].
In [38], the authors calculated the pressure-dependent crystal structures of triclaben-
dazole, an oral anthelmintic. In this work, instead of periodic DFT calculations, the authors
used a “supercell” approach by constructing the 3 × 3 × 3 supercells of two polymorphs.
By comparing the stabilities of two types of TCBZ using the DFT and embedded fragment
techniques and Gibbs free energies comparison, they discovered a phase transition between
Int. J. Mol. Sci. 2023, 24, 14155 22 of 40

TCBZ Forms I and II at 5.5 GPa and ≈350 K. Unfortunately, no such examination of the
pressure-induced phase transition of TCBZ has been previously performed so as to be
compared with the calculated results.

3.3.6. Resorcinol
In the next work [41], the authors have used an unusual approach to describe the
behavior of resorcinol under compression and high temperatures. It was presented that
DFT and DFTB can be coupled in a quasi-harmonic model to correctly present the crystal
structures, phase boundaries, and thermochemistry of the α and β phases of resorcinol.
The authors suggested that the combination of DFT and DFTB+ for phonon calculations
could provide a good balance between the computational time and accuracy of the results
and reduce the calculation time by 1–2 orders of magnitude compared to the more conven-
tional approach based solely on DFT. Although the use of DFTB3-D3(BJ) alone provided
inaccurate results since, at all investigated temperature and pressure values, α-resorcinol
was more thermodynamically stable and no α→β phase transition was observed, the
mixed approach indicates the correct phase behavior in resorcinol. In this case, the results
of Γ point DFT phonon calculations are comparable to the mixed approach. However,
more significant differences could be expected in more complex examples with greater
differences in conformation between polymorphic forms. The revolutionary approach is
based on the phonon density of states calculations in a large supercell with DFTB and
then shifting individual DFTB phonon bands to the frequency obtained from DFT cal-
culations performed in a crystallographic unit cell. As a result, this operation ensures
DFT-quality phonon modes at the Γ point, while the phonon dispersion is modeled with
DFTB. Furthermore, the application of this mixed DFT/DFTB+ method creates the possibil-
ity to compute the thermodynamic properties for larger and more complex systems in a
reasonable time. The QHA B86bPBE-XDM method has been applied to test the pressure
dependence of the lattice constants. In the case of the α phase, the largest mean absolute
error of 1.1% occurs for the b-edge and 1.3% for the a-constant for β-resorcinol, with a good
representation of the experimental data. This model most accurately predicts the lattice
parameters, while, at low pressures, the B86bPBE-XDM and PBE TS methods without
quasi-harmonic thermal expansion underestimate the volumes [54]. Nevertheless, at high
pressure, all approaches overestimate the volumes. From the phase diagram, the authors
estimated the phase transition temperature in the 0–1 GPa regime. At ambient pressure,
the phase transformation was predicted to occur at 368 K, corresponding almost exactly
to the experimental conditions and at 260 K at 0.4 GPa, which is in reasonable agreement
with observations at room temperature. It is noted that the phase transition of resorcinol
strongly depends on temperature and pressure as a kinetic factor plays an important role
in this phase transformation [41].

3.3.7. Glycine
In the next work [42], the authors investigated the structural, electronic, and thermo-
dynamic properties of the α-, β-, and γ-glycine polymorphs. The lattice parameters and the
unit cell volume of the α phase was predicted with a maximum error of 1.42% and 2.53%,
respectively. Based solely on the results of the lattice energy values, the stability order
of the glycine polymorphs was in disagreement with the experimental data. However,
the application of quasi-harmonic calculations allowed them to obtain the correct stability
order. Excellent agreement was found for the isothermal volume behavior of the α-glycine
in the 0 to 50 GPa range. The phase transition of the γ- to α-glycine was recognized at
pressure values of 0.98 and 0.55 GPa at 300 and 400 K, respectively, and the α phase was
indicated to be the most stable at 500 K within the investigated pressure range. The best
predicted value of the γ-to-α-phase transition condition was obtained from a fit to the
third-order Birch−Murnaghan equation of stat and they were 444.55 K and 1 bar, whereas,
experimentally, it was around 440 K. This example has shown that entropy and thermal
contribution could be crucial to accurately foresee the relative stability of the polymorphic
Int. J. Mol. Sci. 2023, 24, 14155 23 of 40

forms and justify the need for Gibbs free energy calculation. In the case of glycine, the
authors also revealed that the calculated order of phases could depend on the dispersion
correction method as the experimental ranking of the lattice energy was obtained by using
MBD dispersion [105] correction, although it was not achieved by using the D3 scheme or
TS correction [106].

4. Fundamental Aspects of DFT Calculations at High Pressure


After the description of some of the most illustrative examples in the previous section,
in this one, a different approach to results presentation is used. Instead of dividing the
reviewed works into particular classes of compounds, below, the fundamental aspects
of DFT calculations are described, starting from structure prediction and optimization
through the calculation of vibrations, which are the essential part of the crystal entropy and
the derivation of the thermodynamic properties. It should be noted here that the solid-state
normal mode calculations are the another level of work compared to the geometry/crystal
cell parameters’ optimization, both in terms of theory and the complexity of calculations
and computational time needed to obtain the results. Finally, the calculations of phonon
properties are described.

4.1. Geomtry Optimization at Various Pressure Conditions and Crystal Structure Prediction
As in other types of studies utilizing molecular modeling methods, in the high-pressure
DFT calculations geometry optimization is a requirement for the further calculation of
various properties. Usually, the SCXRD-determined unit cell is used as the initial basis for
further calculations. However, in the absence of an experimentally determined structure,
the crystal structure prediction methods (CSP) can be used, which are described below.
The goal of the CSP is to determine the crystal structure of a solid based on the
molecular structure. Computational methods that can be used to achieve this goal include
simulated annealing, evolutionary algorithms, distributed multipole analysis, random sam-
pling, basin hopping, data mining DFT, and molecular mechanics calculations. For many
years, the CSP have been investigated as an addition to experimental solid form screen-
ing, assisting in the structural characterization of observed solid forms, understanding
molecular crystallization conditions, and determining when it may be reasonable to stop
screening [107,108]. Moreover, the CSP blind tests—a type of of CSP competition—are very
popular not only among academics, with six editions already successfully completed [109]
and the results of the seventh anticipated to be published [110].
In particular, the CSP methods based on dispersion-corrected density functional theory
(DFT-D), available in the GRACE software, have been successfully used to model multiple
compounds, whose behavior has been then studied also at high pressure [111,112], and this
approach has proven to be very helpful in explaining experimental findings. More such
examples are presented below.
Sometimes, in the case of high-pressure studies, the crystal structure determined
experimentally at ambient conditions is then optimized at higher pressure, which can result
in major structural changes. The presence of structural phase transition can be observed by
abrupt changes in lattice parameters or volume on compression. The authors of the next
work [43] investigated the isosymmetric first-order phase transition of L-histidine. The DFT
calculations were performed with symmetry constraints. As a result, the discontinuities in
the lattice parameters and unit cell volume showed the formation of new phases I0 and II0 .
The theoretical phase transformation for the orthorhombic form (I) was observed between
4.45 and 4.62 GPa, which is in good agreement with the experimental transition pressure
equal to 4.5 GPa. In the case of the monoclinic phase, a reduction in volume was observed
ca. between 3.00 and 3.21 GPa and experimentally at 3.1 GPa.
Discontinuities in lattice parameters were also observed in the next work [65] where
the authors investigated the bis-1,2,3-thiaselenazolyl radical dimer under compression up
to 14 GPa. Surprisingly, the isostructural radical dimer, where the hydrogen atom replaces
fluorine, has remarkably different behavior under higher pressure despite the similarity
Int. J. Mol. Sci. 2023, 24, 14155 24 of 40

in the structure and transport properties at ambient pressure. The dimer with hydrogen
[1a]2 undergoes phase transformation at around 5 GPa, with good agreement with the
experimental part of the work, whereas the fluorine dimer undergoes uniform compression
without phase transition. The structural evolution and related changes in the molecular
and band electronic structures of the two compounds were impressively replicated by DFT
calculations. In the case of [1a]2 , abrupt changes were observed in the lattice parameters
and bandgap values as an effect of pressure. From enthalpy difference calculations, it was
projected that, at low pressures, the uniformly compressed σ-dimer of the [1a]2 structure
is more stable, with the buckled π-dimer taking dominance at higher pressures. The
pressure transition was calculated to be in good accordance with the experimental results.
The authors run two sets of calculations of [1a]2 dimers as geometry optimization was
performed at 0 K and the spontaneous phase transition from the σ- to π-dimer was not
expected due to the energy barrier. To overcome this issue, for the π-dimer form, the
observed crystal geometry at different pressure values was used as an initial structure for
calculations. The second set included isostructural compression of the ambient pressure
crystal geometry of the ambient structure of [1a]2 to generate the equation of state for the
hypothetical σ-dimer variant of [1a]2 . For [1b]2 , the crystal geometry at various pressure
was applied as a starting point for geometry optimizations.
In the last work of this section [47], the aim of the study was to investigate the high-
pressure polymorphism of L-threonine within the pressure ranging from ambient pressure
to 22 GPa using single-crystal X-ray and neutron powder diffraction. Due to the limited
scattering geometry of the diamond anvil cell, the high-pressure study suffered from
low completeness. To overcome this issue, the bond distance and angles were restrained
to values observed at ambient pressure. The periodic DFT calculations were applied to
validate the high-pressure structural models. Therefore, the geometry optimizations were
performed on single-crystal structures at ambient pressure and those up to 22.31 GPa with
cell dimensions and space group symmetry restraints. From a comparison of the bond
distance and angle changes, it was detected that the torsion angles did vary significantly
with pressure. Experimentally, it was noticed that the compound exhibited three phase
transitions triggered by increasing pressure. DFT calculations were also used to calculate
the pressure dependence of the lattice energy of L-threonine. The lattice energy increased
gradually in phases I, I0 , and II, and a discontinuity occurred during the II to III phase
transition at 18.2 GPa, resulting in the abrupt destabilization of the lattice energy.
As stated above, CSP methods have found application also in DFT calculations at high
pressure. Currently, there are multiple available codes that enable such predictions, i.e.,
GRACE, Polymorph Predictor, or USPEX. This last one has gained particular attention
in works describing DFT calculations at high pressure. Universal Structure Predictor:
Evolutionary Xtallography (USPEX) is an algorithm and software package developed
for the prediction of crystal structures. It employs an evolutionary line of action in the
investigation and prediction of the stable or metastable crystal structures of compounds.
The algorithm’s ability to propose novel crystal structures has proven to be a great asset for
material discovery and understanding structure–property relationships, and a great help in
guiding experimental synthesis works.
In the first work [23], the authors used the USPEX algorithm to find the candidate
high-pressure structure of ethylenediamine bisborane (EDAB); see Figure 6. As a result,
they obtained 825 structures and, to shorten the list, they also compared the experimental
and simulated XRD patterns of these structures. The 14 most promising structures for phase
II (at 1.72 GPa) and 17 candidate structures for phase III (at 11.44 GPa) were subsequently
optimized using the nonlocal van der Waals density functional theory (vdW-DF). This is
regarded to give more accurate relative enthalpies than the semilocal GGA PBE functional.
Following this, the calculated differences in enthalpies and parameters allowed the authors
to establish the final top seven structures. Simultaneously, they performed calculations
with the use of the PBE functional, which led to the same sequence of the lowest-enthalpy
structures in the 0−17 GPa range, although with different values of transition pressure.
Int. J. Mol. Sci. 2023, 24, 14155 25 of 40

However, DFT calculations suggest that phase II is P21 /c (A) and phase III is P21 /c (B),
irrespective of the functional used. To confirm these identifications, Rietveld’s refinement
method was employed to enhance the experimental XRD patterns for phases I−III using
the corresponding theoretical crystal structures. Based on a comparison of the experimental
and theoretical structures, the agreement for phase I was almost perfect. As compressed
EDAB could be a potential candidate for hydrogen storage, its mechanical properties
have been investigated at high pressure. The discontinuity of the cell parameters on the
pressure reveals two phase transitions at around 1 and 7 GPa; however, the second drop
is relatively small, suggesting the coexistence of phases II and III, which is supported by
experimental results. In addition, the authors performed phonon calculations using solely
the PBE functional for the PBE-optimized parameters since vdW-DF was not accessible
in the version of Quantum Espresso used. The EDAB structures classified as phases I−III
Int. J. Mol. Sci. 2023, 24, x FOR PEER REVIEW 26 of 41
were confirmed to be true minima by the phonon calculation as no negative frequencies
were detected.

Figure6.6.Pressure-induced
Figure Pressure-induced phase
phase transitions
transitions of crystalline
of crystalline ethylenediamine
ethylenediamine bisborane bisborane [23].
[23]. Copyright
Copyright 2023 American Chemical
2023 American Chemical Society. Society.

In the
In the next
next work
work [59], [59], the
the authors
authors used
used USPEX
USPEX algorithms
algorithms to to generate
generate candidate
candidate
structures of
structures ofcroconic
croconicand andsquaric
squaric acid
acid under
under 25 GPa.
25 GPa. In order
In order to establish
to establish the phase
the phase tran-
transition, the pressure dependence of the O-H bond
sition, the pressure dependence of the O-H bond under compression and decompression under compression and
decompression
was examined using was examined using DFT
DFT calculations. In calculations. In the case
the case of croconic andof croconic
squaric and
acid, squaric
isotropic
acid, isotropicpredicts
compression compression predicts of
the formation thethe
formation of the
crystalline crystalline
structure at 12structure
GPa andat15 12GPa,
GPa
and 15 GPa,The
respectively. respectively.
formation The of a formation
new polymeric of a structure
new polymeric structure
of croconic of croconic
acid occurs acid
at around
occurs
15 at around
GPa under 15 GPa
isotropic under isotropic
compression. In addition,compression.
the outcomes In of
addition, the outcomes
evolutionary simulations of
evolutionary
at 25 GPa revealed simulations at 25 GPaofrevealed
that the structures both acidsthat the structures
are comparable of both acids
to structures obtained are
under isotropic
comparable compression
to structures in terms
obtained of enthalpy.
under isotropicUpon decompression
compression in terms forofboth acids,
enthalpy.
higher pressure stabilizes
Upon decompression for the
both structure to lower
acids, higher pressure.
pressure It was noticed
stabilizes that the
the structure to phase
lower
transition
pressure. It of was
squaric
noticedacidthat
could thebephase
foundtransition
under isotropic compression
of squaric acid couldand USPEXunder
be found evo-
lutionary simulations and
isotropic compression anditUSPEX
was associated with simulations
evolutionary the formation andof O-H-O weak covalent
it was associated with
bonds, and theofhigh-pressure
the formation O-H-O weak form covalenthadbonds,
the highandsymmetry group I4form
the high-pressure m. Based
had theonhigh
the
calculated dispersion phonon curves, the lack of negative (imaginary)
symmetry group I4 m. Based on the calculated dispersion phonon curves, the lack of frequencies indicates
that structures
negative obtained
(imaginary) at 25 GPaindicates
frequencies are stable forstructures
that both croconic and at
obtained squaric
25 GPaacids. The
are stable
authors also compared the simulated Raman spectra with the
for both croconic and squaric acids. The authors also compared the simulated Raman experimental data.
Another
spectra with example where thedata.
the experimental evolutionary structural search algorithm USPEX was applied
to predict
Anotherstableexample
polymorphic where structures at high pressure
the evolutionary was described
structural by ClarkeUSPEX
search algorithm et al. inwas
the
case of α-glycylglycine (α-digly) [46]. From the results of their calculations,
applied to predict stable polymorphic structures at high pressure was described by Clarke the orthorhombic
phase
et al. inisthe
suggested to be the lowest-enthalpy
case of α-glycylglycine (α-digly) [46]. polymorph above 6.4
From the results GPa.calculations,
of their However, the the
discrepancy between the measured and simulated PXRD
orthorhombic phase is suggested to be the lowest-enthalpy polymorph above implies that the orthorhombic
6.4 GPa.
structure
However, (digly-P2 1 21 21 ) wasbetween
the discrepancy not experimentally
the measured observed
and and the monoclinic
simulated phase persisted
PXRD implies that the
orthorhombic structure (digly-P212121) was not experimentally observed and the
monoclinic phase persisted for the full pressure range. Both Raman and PXRD data
indicated that the P21/c symmetry was maintained, suggesting an isosymmetric phase
transition at around 6–7.5 GPa. Hence, the transition to the orthorhombic phase requires
a change in packing and bending peptide backbone; the transformation to this form might
Int. J. Mol. Sci. 2023, 24, 14155 26 of 40

for the full pressure range. Both Raman and PXRD data indicated that the P21 /c symmetry
was maintained, suggesting an isosymmetric phase transition at around 6–7.5 GPa. Hence,
the transition to the orthorhombic phase requires a change in packing and bending peptide
backbone; the transformation to this form might not have been noticed in the experiment
due to the large energy barrier. An isosymmetric phase transition from α-digly to α0 -digly
was predicted at 10 GPa, as ZPE was not included in the enthalpy, or 11.4 GPa, when the
ZPE was added, whereas the experimental transition pressure was ca. 6.7 GPa. However, it is
worth mentioning that the differences in enthalpy in the 6–9 GPa range were smaller than
10 meV/atom. In this work, the thermal contribution to the energy was not included regarding
the fact the supercell approach should be applied, which is very computationally expensive.
Adding the ZPE is supposed to give results in better agreement with the experiment, with
a simultaneously lower computational cost. Nevertheless, the changes in the a- and c-axes’
lengths under compression agreed with the experimental observations. The largest differences
were observed for the b-axes. This could be associated with the fact that the b-direction
includes multiple H-bonding interactions, which are known to be difficult to accurately model
using the PBE functional.
The USPEX algorithm was also used for the prediction of structures of 2,4,6-trinitro-
3-bromoanisole (TNBA), an energetic material described in the previous section, “High-
Energetic Organic Materials”.

4.2. Vibrational Spectra


Pressure-induced phase transition can not only affect structure parameters but also
modify the vibrational spectra obtained for studied compounds. Therefore, the examina-
tion of vibrational modes by in silico methods could potentially bring great insights into
the analysis and assignment of spectral features to different molecular motions, such as
stretching, bending, or torsional vibrations, as the pressure changes.
Confirmation of this statement can be found in the work about glycinium maleate [40];
the DFT calculations can be used to assign inter- and intramolecular modes and help to
describe local changes in molecules under compression. Based on theoretical research, the
band evolution of glycinium maleate associated with the internal modes strongly indicated
that the conformation was more correlated with the maleic acid molecules than with the
glycine molecules. The DFT calculations were particularly useful in the assignment of the
bands between 3050 and 3100 cm−1 due to extensive changes in the intensity of the bands
under compression.
The authors of the next work [58] investigated the behavior of chloroform and
dichloromethane upon compression and suggested the type of interaction that governed
the compressibility. They found that for both compounds, the most stable phase was the
ambient-pressure structures (Pnma and Pbcn for CHCl3 and CH2 Cl2 , respectively), which
persisted up to 32 GPa, and the calculated phonon dispersion curves did not exhibit imagi-
nary frequencies, suggesting the stability of these structures at high pressure. Furthermore,
their results indicate that the P63 structure of CHCl3 is metastable. Comparison of the calcu-
lated results with the experimental information of the Raman band positions allowed them
to characterize the pressure-induced evolution of the crystal structures of both compounds.
In the next work [60], DFT calculations were used to simulate the Raman spectra of
a single molecule and crystal of diisopropylammonium perchlorate (DIPAP) to assign and
compare the calculated results with experimental measurements. In the case of molecular
calculations, it was noticed that including the dispersion correction had no significant impact
on the vibrational frequencies in the region 200−1650 cm−1 . However, there was a slight
shift in the higher wavenumber mode. For solid-state calculations, the results were in good
agreement with experimentally obtained intermolecular and intramolecular vibrational modes
as dispersion correction has a minor impact on the accuracy of theoretical results.
In the next study [63], the authors investigated the pressure dependence of the Raman
spectra of sorbic acid crystals to acquire insights into their vibrational properties. The
observed broadening of the Raman bands due to the decrease in the cell volumes and
Int. J. Mol. Sci. 2023, 24, 14155 27 of 40

higher disorder, the upshift of the wavenumbers, decreases in the Raman intensities, and
some discontinuities were the main observed differences as an effect of the applied pressure.
It is worth mentioning that at 0–6.0 GPa, which is in the range of expected transition
pressure, the pressure dependence wavenumbers were nonlinear. The alterations detected
in the Raman spectra’s external mode area were interpreted as changes in hydrogen
bonding, resulting in changes in the conformation of the molecules. This interpretation
was supported by the observation that the main changes observed among the internal
modes were associated with the bands of a carboxylic group that actively took part in
the hydrogen bonds. Based on the results, it was stated that the sorbic acid underwent a
conformational transition due to the similarity of the Raman spectra between the room-
temperature monoclinic structure and the high-pressure phase.
In the next work [54], the authors investigated the α, β phases of resorcinol by combining
Raman, time-domain terahertz, and inelastic neutron scattering spectroscopy with solid-state
density functional theory (DFT) calculations. The crystal structures of both polymorphs were
optimized by using various GGA functionals to choose the most accurate by calculating the
mean absolute errors between experimental and calculated internal coordinates. Based on the
results, PBE TS was chosen for the examination of the lattice constants of both phases under
compression. Even at high pressure, the deviations from the predicted values were less than
3%, except for the b-constant in the β phase, which was observed to decrease more easily.
According to the vibrational spectra analysis, it was indicated that the standard (PBE) and
“hard” (rPBE) GGA functionals were able to successfully reproduce the whole spectral range,
aiding in the detailed interpretation of the bands.
In the next study [36], the authors investigated the structural and vibrational properties
of 1,3,5,7-cyclooctatetraene (COT) under pressure. The calculated Raman-active modes of a
single COT molecule were in good agreement with experimental values. The theoretical
lattice constants obtained from ab initio molecular dynamics simulation at 0 GPa and
3.8 GPa were in excellent accordance with measurements.

4.3. Enthalpy (∆H) Calculations


In most cases, the transition pressure has been calculated based on differences in
enthalpy between polymorphs, since the calculations are conducted at 0 K, where the
enthalpy is equal to the Gibbs free energy. The phase transformation occurs at a pressure
value where ∆H is equal to zero. As shown in the examples below, the predicted transition
pressure from ∆H calculations can agree well with experimental findings [48,50]. However,
there have been reported cases where a lack of entropy and thermal contribution could lead
to the incorrect order of stability of the polymorphic phases, as in the case of glycine [94], or
not even indicate the experimentally observed phase transition, as for the α and β phases
of resorcinol [41].
In the first study [48], the authors performed dispersion-corrected DFT calculations of
the α and β phases of L-threonine to investigate their structural, electronic, and dielectronic
properties at hydrostatic compression. The root mean square deviation (RMSD) of the lattice
constant between the calculated and experimental values was ca. 1% and the difference in
intra- and intermolecular bond lengths and angles was around 6–8%, corresponding to the
experimental results. From the differences in enthalpy (∆H) values, the phase transition
from the α and β phases was observed at 1.5 GPa, which agreed with the experimental work.
On the pressure dependence of lattice parameters, a discontinuity at 1.5 GPa was noticed,
suggesting the occurrence of a phase transformation. The analysis of changes in the lattice
parameters of intra- and intermolecular bonds revealed that the changes in lattice structure
played a predominant role in the pressure effect on intermolecular hydrogen bonds, but
not the change in intramolecular bond lengths. Furthermore, it was established that the β
form could exist without pressure as a metastable phase. In the next work [57], Howard
et al. studied the phase transition and thermal expansion of ammonium carbamate. This
compound is known to exist in two polymorphic forms, namely α and β. By employing
DFT calculations, the authors showed that the thermodynamically stable phase at ambient
Int. J. Mol. Sci. 2023, 24, 14155 28 of 40

pressure is the α-polymorph, with a calculated enthalpy difference with respect to the
β-polymorph of 0.399 kJ mol−1 . According to the calculation results, the transition to the
β-polymorph could occur at 0.4 GPa. However, this observation has not yet been studied
using experimental methods. Another interesting aspect of this work was that the authors
considered both the dispersion-corrected (TS) and non-corrected DFT calculations. For
the work reported, they found that geometry optimizations completed without dispersion
corrections led to disagreements of 10% in the unit cell parameters, these being reduced
to 1–3% with the use of dispersion corrections. Similarly, the estimated phase transition
pressure reported later was reduced from 2.0 to 0.4 GPa. In contrast to other similar works
reported, when determining the pressure of phase transition, the authors compared solely
the ∆H values instead of ∆G. This was due to the calculations being performed solely at
0 K, without taking into account either ZPVE or entropy changes.
In another work [55], the authors used DFT methods to investigate the structural
changes of glycine and L-alanine crystals under pressure up to 10 GPa. They performed
calculations with and without van der Waals interaction correction. Based on the results,
the authors established that the inclusion of van der Waals interactions is crucial for an
accurate description of the intermolecular interactions of amino acids. The differences
in the stability of glycine polymorphs are very small, so it is challenging to accurately
establish the stability order. Furthermore, the calculated differences are sensitive to the
exchange–correlation functionals applied. In the case of glycine, the calculations with the
vDW interaction better predicted the stability order of polymorphs, unlike calculations
without this component. Based on the calculated enthalpy values, the pressure-induced
phase transition of the β to δ phase was predicted at 4.0 GPa, 1.5 GPa, and 0.7 GPa for PBE,
vdW-DF, and vdW-DF-c09x, respectively, whereas, experimentally, it occurs at 0.76 GPa.
Consequently, it was stated that PBE without including dispersion correction overestimated
the transition point. For L-alanine, only the vDW-aware functionals correctly determined
the unit cell parameters under pressure.
In the next work [56], the DFT and PIXEL methods were used to research the behavior
of L-serine at pressure up to 8.1 GPa. The experimental structure parameters were reopti-
mized with fixed unit cell dimensions and unrestricted other parameters and symmetry.
The obtained geometries were in good agreement with the neutron powder diffraction
analysis for each phase, and the root mean square deviation (RMSD) of the bond length
in the structures was never greater than 0.08 Å. Based on the results, performing the ge-
ometry relaxation of an ambient-pressure form with a fixed external pressure parameter
could be a powerful tool for the prediction of the high-pressure-induced differences in
molecular packing and geometry. The PIXEL method is based on the determination of a
molecular electron density map, followed by processing the map into larger pixels and
then the calculation of energy terms between pairs of pixels in adjacent molecules. The
technique does not take into account any intramolecular changes in energy and it is used to
calculate the energies of intermolecular interactions. However, the transition from phase
I to phase II is partly driven by a conformational change in the molecule. According to
PIXEL calculations, there is a significant energy gap between the intermolecular energies
of phases I and II. The phase transition from L-serine I to L-serine II occurs near 5 GPa,
where there is a break in the gradient of the curve, with the enthalpy of phase II being more
negative after the change than the extrapolated values for phase I due to the stabilization
of the internal energy of the serine molecules and volume reduction. It seems that the II to
III phase transformation is driven by rearrangement in intermolecular interactions.
In a subsequent paper examining the behavior of L-serine crystals at high pressure [50],
the authors applied DFT calculations to investigate the mechanism of reversible phase
transition between polymorphs of small organic molecules like L-serine. Firstly, they per-
formed optimizations of the crystal structures documented to dominate at these pressures.
However, this approach did not allow the researchers to compare the enthalpies and crystal
energies of different phases at the same pressure and explain the occurrence of phase
transition. Hence, the authors performed optimization of all polymorphs at the same
Int. J. Mol. Sci. 2023, 24, 14155 29 of 40

pressure range, outside of the range of their existence. Based on the enthalpy differences
(∆H) between polymorphs, the pressure transitions were defined as 3.7 GPa and 5.3 GPa
for the transition from phase I to phase II and phase II to phase III, respectively. Although
Int. J. Mol. Sci. 2023, 24, x FOR PEER REVIEW 30 of 41
the calculated values were smaller than the experimentally registered 5.3 GPa and 7.8 GPa,
respectively, the difference between the points of the two-phase transitions was close to
that observed
a wide experimentally.
pressure range; The discrepancy
thus, the transitions could be might occur
initiated at due to thepressure
different limitations
values,of the
depending on the experimental conditions. Although the entropy factor was not included, of
model or might suggest kinetic control over the transition. Additionally, the polymorphs
theL-serine have been
computational experimentally
simulations observed
correctly predictedto coexist in aofwide
the order pressure range;
the polymorphs’ thus, the
stability.
transitions could be initiated at different pressure values, depending on the experimental
As, during phase transitions I–III, the crystal symmetry did not change significantly and
conditions. Although the entropy factor was not included, the computational simulations
no crucial changes were observed in the Raman spectra, the authors assumed the
correctly predicted the order of the polymorphs’ stability. As, during phase transitions I–III,
differences in the TΔS term into Gibbs free energy to be small and not including this factor
the crystal symmetry did not change significantly and no crucial changes were observed in
allowed them to reduce the computational time, with minor inaccuracy for the Gibbs
the Raman spectra, the authors assumed the differences in the T∆S term into Gibbs free
energy. Based on the research, the transition of L-serine polymorphs was shown to be
energy to be small and not including this factor allowed them to reduce the computational
triggered by the PV term, at least between the I and II phase. The investigation of the effect
time, with minor inaccuracy for the Gibbs energy. Based on the research, the transition of
of pressure on H-bond interactions in gas-phase cluster models helped to highlight the
L-serine polymorphs was shown to be triggered by the PV term, at least between the I and
significant difference between the I→II and II→III phase transitions. The first transition is
II phase. The investigation of the effect of pressure on H-bond interactions in gas-phase
triggered by the large overstrain of a selected intermolecular hydrogen bond,
cluster models helped to highlight the significant difference between the I→II and II→III
experimentally manifested in changes in cell parameters, and the large hysteresis; the
phase transitions. The first transition is triggered by the large overstrain of a selected
latter is proposedhydrogen
intermolecular to be accompanied by multiple
bond, experimentally small changes
manifested in various
in changes in cellhydrogen
parameters,
bonds. The current study demonstrates that relatively simple
and the large hysteresis; the latter is proposed to be accompanied by multiple small calculations withchanges
the
combination of detailed experimental data may provide insight
in various hydrogen bonds. The current study demonstrates that relatively simple calcu-into the macro- and
micro-driving
lations withforces of pressure-induced
the combination of detailedphase transition data
experimental in hydrogen-bonded
may provide insight molecular
into the
crystals.
macro- and micro-driving forces of pressure-induced phase transition in hydrogen-bonded
molecular crystals.
4.4. Gibbs Free Energy (ΔG) Calculations
4.4.InGibbs
orderFreetoEnergy
obtain (∆G)
theCalculations
thermodynamic properties, the very computationally
demanding andto time-consuming
In order phonon properties,
obtain the thermodynamic density ofthe states calculations need
very computationally to be
demanding
performed. Therefore, phonon
and time-consuming the size density
of the of
investigated crystalsneed
states calculations and to the
berequired
performed. computer
Therefore,
parameters
the size ofcould limit the application
the investigated of this
crystals and the approach to determine
required computer the relative
parameters couldstability
limit the
and transitionofpressure.
application Furthermore,
this approach to determinethe calculation of thermodynamic
the relative stability and transition functions
pressure. canFur-
thermore,
help the calculation
to establish whether the of thermodynamic
phase transitionfunctions can help
is enthalpy- (e.g.,toglycine
establish[49])
whether the phase
or entropy-
transition
driven (e.g., is enthalpy- (e.g.,[39],
chlorothiazide glycine [49])7).
Figure or entropy-driven (e.g., chlorothiazide [39], Figure 7).

Figure 7. Thermodynamic calculations for chlorothiazide, indicating entropy-driven phase transition


Figure 7. Thermodynamic calculations for chlorothiazide, indicating entropy-driven phase
of FormofI to
transition Form
Form I toIIForm
at higher
II at pressure. Adapted
higher pressure. from [39],
Adapted fromlicensed under CC
[39], licensed BY 4.0.
under CC BY 4.0.

The study of glycine [49] is one example where the Gibbs free energy (ΔG) was
calculated to predict phase transformation. The aim of the work was to evaluate whether
periodic DFT calculations could be applied to the investigation of the dependence of
increasing pressure on the molecular crystal structures of glycine and their stability.
Firstly, the authors tested various functionals with or without dispersion correction to
choose the most accurate based on the results of the geometry optimization of the γ and δ
Int. J. Mol. Sci. 2023, 24, 14155 30 of 40

The study of glycine [49] is one example where the Gibbs free energy (∆G) was
calculated to predict phase transformation. The aim of the work was to evaluate whether
periodic DFT calculations could be applied to the investigation of the dependence of
increasing pressure on the molecular crystal structures of glycine and their stability. Firstly,
the authors tested various functionals with or without dispersion correction to choose the
most accurate based on the results of the geometry optimization of the γ and δ phases
at 3.27 GPa, as, at this pressure value, both forms have been observed experimentally.
According to the outcomes, the GGA functionals provided more accurate results than
LDA. Surprisingly, including OBS dispersion correction decreased the accuracy of the LDA
calculations. For further investigation, the PBESOL functional was chosen, since it was
found to be very accurate in most cases, but also was originally developed for densely
packed solids. The next step was to perform geometry optimization up to 7.8 GPa for
both polymorphs. The obtained results of the lattice parameters were found to be in very
good agreement with the corresponding experimental data. To test the accuracy of the
calculations, the PBE functional without dispersion correction at pressure values of 0 GPa,
5.83 GPa, and 7.80 GPa was applied. From the comparison of the results obtained using
two different functionals, the crystal structures were significantly less accurately modeled
using PBE than PBESOL. Additionally, the authors suggested that the simulated lattice
parameters in some cases could be even more accurate than experimentally measured as
they better fit the overall trend of changes in lattice parameters within increasing pressure.
Performing thermodynamic parameter calculations of the γ and δ phases at various values
of pressure allowed them to establish the order of stability and the transition pressure at
which the order was reversed. The calculated values of Gibbs free energy (∆G) were in
good accordance with the experimental findings, where the γ phase started to transform
into the δ polymorph at ca. 2.74 GPa. Furthermore, based on the computational results,
the phase transition of glycine was described as enthalpy-driven. Nonetheless, the entropy
effect was also favorable for the transformation to occur at high pressure.

4.5. Phonon Calculations


Since, in some cases, the predicted high-pressure structures could be metastable
phases, to test their stability, the phonon frequency can be calculated to detect the presence
of potential imaginary frequencies [23,58,59,61]. The lack of negative frequencies suggests
that the forms predicted by theoretical methods are stable under the investigated conditions.
For instance, no imaginary frequencies were observed for the structures of croconic and
squaric acid obtained at 25 GPa for both acids [59], or in the case of chloroform and
dichloromethane for the structures obtained at 32 GPa. However, the previously proposed
structure P63 of CHCl3 was found to be a metastable form; thus, the ambient-pressure
structure is likely to be a ground-state polymorph of this compound up to 32 GPa [58].
Comparably, based on phonon calculations, the structures of EDAB identified as the I–III
phases were proven to be true minima [23]. The phonon dispersion curves were evaluated
along the high-symmetry direction to obtain information about the structural stability
of the g, α, d, and β phases of C11 N4 . According to the results, the g phase was the
most energetically stable at ambient pressure, and α- and d-C11 N4 indicated no negative
frequencies, in contrast to β- and g-C11 N4 [61].

5. Other Aspects Associated with DFT Calculations at High Pressure


In this section, the problems and phenomena typical of DFT studies at high pressure are
presented. These aspects include methods of determination of the pressure-induced phase tran-
sition conditions, energy barrier value calculations, situations in which no pressure-induced
transformations are observed, or applications of computational anisotropic compression.

5.1. Determination of Pressure-Induced Phase Transition Conditions


There are several known methods to determine the transition pressure based on
structural, energetic, or property changes. Below, some of them are listed with a short
Int. J. Mol. Sci. 2023, 24, 14155 31 of 40

description of the work, to give an overview of the methods used and where they have
been applied. Notably, the mentioned approaches are not limited to the given examples,
and some of them, solely or in combination with different approaches, have been used in
other papers.

5.1.1. Common Tangent to the Two E(V) Curves, p = −dE/dV


As shown in the case of C11 N4 study [61], the transition pressure can be obtained
directly from the calculated energy–volume curves. The formation enthalpy under the
pressure is defined as H = E (V) + PV, where E is the total energy for the cell with volume
(V). The transition pressure for each phase transformation was obtained by calculating the
common tangent slope of the two energy–volume curves. It was shown that the transition
from g-C11 N4 to α-C11 N4 , d-C11 N4 , and β-C11 N4 can occur at 3.557 GPa, 9.468 GPa, and
46.032 GPa, respectively. The obtained data indicated that the transition from g-C11 N4
to α-C11 N4 and d-C11 N4 occurs at low pressure, and the transition from g-C11 N4 to β-
C11 N4 occurs in the high-pressure phase. Additionally, the authors investigated the phase
transitions and vibrational, mechanical, and thermodynamic properties of four polymorphs.
Thanks to the performed calculations, the prediction of the mechanical and thermodynamic
properties, including the bulk modulus, heat capacity, and thermal expansion coefficient,
of C11 N4 polymorphs was possible.

5.1.2. Changes in Properties Observed upon Compression


In some papers, the possibility of pressure-induced phase transition was suggested
based on changes in properties, e.g., upon noticing alterations in bandgap widths for ε-CL-
20, 400 GPa was supposed to be the critical pressure for insulator–metal transition [37], or
in the case of the previously described LLM-105, the abrupt bandgap decrease suggested
a pressure-induced phase transformation, confirmed by experimental work [24]. The
electronic properties can change with an increase in pressure, which can be used in tuning
materials’ properties by exposing them to a high pressure factor [113]. Consequently, the
theoretical investigation of electronic and mechanical properties can be a beneficial tool in
designing materials with desired properties, e.g., for photovoltaic applications.
In one work [64], single-point DFT band structure calculations were performed to
investigate the pressure-induced changes in the electronic structures of the α, β, and γ
phases of the oxobenzene-bridged 1,2,3-bisdithiazolyl radical conductor (3a) obtained
at 0 GPa (α phase), 6.0 GPa (β phase), and 11.1 GPa (γ phase). The α phase is a Mott
insulator, as the bandwidth of all explored bands was small. The calculated results showed
that the increase in pressure resulted in an increase in bandwidth in the β phase; further
compression initiated the second phase transition to the γ phase, lifting the HOMO band
and causing higher-lying virtual orbitals to drop. Generally, the presence of low-lying
LUMO triggers high electron affinity and creates an electronically much softer radical
with a low on-site Coulomb potential U, providing an important insight into the design of
radical-based conductors.
In the study of the polymer [Zn(µ-Cl)2 (3,5-dichloropyridine)2 ]n [13], as was mentioned
before, the transition associated with breaking symmetry from P4b2 to P4 was predicted
by using DFT calculations based on the pressure-induced softening of low-frequency
vibrations in the Raman spectra [13]. Therefore, calculations can support the interpreta-
tion of experimental spectra and help to detect structural phase transitions, although the
pressure–volume (p-V) curves may not show any discontinuities.

5.2. Lack of Pressure-Induced Phase Transition


Although phase transition has not been observed in all studies, DFT calculations
have allowed us to simulate the effect of high pressure on the structural, mechanical, and
electronic properties of investigated materials. The lack of structural phase transition was
observed, e.g., in the study of TATB [28], LLM-105 [31,32], and silver fulminate [34] or
aspirin and paracetamol, as described below [53].
Int. J. Mol. Sci. 2023, 24, 14155 32 of 40

In one study [53], the authors investigated aspirin and paracetamol’s (Forms I and
II) behavior up to 5 GPa using long-range dispersion-corrected hybrid density functional
calculations. In both cases, the phase transition was not observed, which was in line with
the experimental results of paracetamol, and no investigational results for aspirin were
available to compare at that moment. In another investigation, it was noted that the phase
transition between Form I and Form II of aspirin was not detected up to 10 GPa using
micro-Raman spectroscopy but, instead, Form I was transformed into a new phase, Form
III, above ~2 GPa, and the authors proposed Form III to be the most stable polymorph of
aspirin at high pressure.
Based on the pressure dependence of the lattice parameters and bond length of both
polymorphs of paracetamol, the maximum compression was observed in the b-lattice
direction, perpendicular to the hydrogen-bonded molecular layers. In the case of Form
I, unusual behavior was observed as compression resulted in the expansion of the a-
lattice parameter. Firstly, the a-edge decreased from 12.621 Å to 12.471 Å at 1.5 GPa
and then increased up to 12.738 Å at 4.0 GPa. Although this is not ordinary behavior
under compression, the expansion of the a-lattice parameter of Form I was experimentally
observed at a pressure range of 2–4 GPa.
From the research on aspirin, both forms were characterized by similar anisotropic
strain at applied pressure. Furthermore, the IR spectra was calculated for each form at
0 GPa and 5 GPa as changes in intermolecular interactions manifested in the shifting of the
bands. The pressure-induced structural alterations were complemented by bands with red
shifts in the IR spectra of the four investigated forms.
The behavior of aspirin crystals under pressure was also an area of interest in a
previously mentioned work [52], where the energy barrier was calculated, indicating that
no pressure-induced phase transition was observed.

5.3. Anisotropic Compression


Anisotropic compression can be simulated by using DFT calculations. It is especially
beneficial in the investigation of energetic materials, as detonation is a non-equilibrium,
ultrafast process with strong orientation dependence, and real-time measurements of this
type of material are challenging. The pressure-dependent vibrational frequency shifts
and the possibility of phase transformation under various compression orientations could
remain elusive [25]. Therefore, a theoretical approach can support the examination and
give valuable insights into the behavior of materials under non-hydrostatic pressure, as in
the study of energetic materials FOX-7 [25], β-HMX [33], and MAPbBr3 [17] or polymer
Zn(µ-Cl)2 (3,5-dichloropyridine)2 ]n [13] and as described below for oxalic acid [62].
In their work [62], the authors analyzed dihydrate and the α and β polymorphic
forms of anhydrous oxalic acid. The calculated structural parameters were in good agree-
ment with the experimental data. The sets of geometry optimization were performed
at several different pressure values applied along the direction of the minimum Poisson
ratio. According to the obtained results, it was shown that the dihydrated form of oxalic
acid undergoes a phase transition as an effect of negative pressure induction smaller than
approximately −0.045 GPa. Both forms of anhydrous oxalic acid also undergo pressure-
induced phase transitions, in this case for positive pressures, larger than around 1.91 GPa
for the α polymorph and around 0.21 GPa for the β polymorph. It was noticed that α- oxalic
acid underwent a pressure-induced phase transition under the effect of applied pressure
directed along the (0.01, 0.73, 0.68) direction, whereas for the β phase, this occurred along
the (0.95, 0.00, 0.31) direction.

5.4. Polymorphic Transition Energy Barrier Calculations


In the case of indole, as well as other compounds, it has been shown that, sometimes,
geometry optimization, even at a higher level of theory, is not sufficient to observe phase
transition. In one work [66], the authors investigated indole crystals’ responses to hydro-
static pressure using molecular dynamics and DFT calculations. Both methods effectively
Int. J. Mol. Sci. 2023, 24, 14155 33 of 40

reproduced the experimental structure with an error of around 1%. However, abrupt
changes in lattice parameters were not observed in the DFT calculation results, even though
the HB to β phase transition energy was relatively small (around 0.1 eV per unit cell).
Nevertheless, the transition barrier was ca. 0.9 eV per unit cell and required increasing
kinetic energy, which was not provided in the static calculations performed at 0 K. These
results revealed that the pressurization was not sufficient for the transformation to occur.

6. Conclusions
Density functional theory is a potent and often used quantum mechanical tool for the
examination of different features of matter. The studies in this area cover a wide variety
of topics, including the creation of original analytical methods centered on the creation of
exact exchange-correlation functionals and the application of this method to the prediction
of the molecular and electronic configurations of atoms, molecules, and solids in both gas
and solution phases. Since there are still problems to be solved, designing and evolving
more effective density functionals is a continual process. Ensuring the appropriate qualities
at a reasonable processing cost is a major quantum task. Future work will concentrate on
creating even more consistently accurate density functionals for particular applications,
enabling researchers to utilize DFT’s relatively high precision at a low cost.
The creation of precise force fields generated from first-principles data is perhaps one
of the most significant directions in the first-principles modeling of molecular crystals.
The types of molecular dynamics simulations and multiscale models required to properly
comprehend the thermodynamics and kinetics of molecular crystals as a function of temper-
ature and other factors would be made possible by such force fields. A complete picture of
the creation, stability, and characteristics of molecular crystals may be obtained by properly
determining the underlying electronic energies and significant response qualities using
entirely first-principles approaches.
The above review shows that DFT calculations can be successfully used to describe
the various molecular solids, including organometallic compounds and high-energetic
materials, which are particularly demanding and hazardous from the experimentalist’s
point of view. Moreover, the types of calculations and modeled properties are diverse and
include not only the most popular geometry optimization but also the ab initio molecular
dynamics, spectral (Raman, UV–Vis, IR), or thermodynamic property calculations. There-
fore, the application of DFT methods can lead to a better understating of the changes in the
structures and properties of materials resulting from the application of high pressure.
Various strategies are being used to determine whether the pressure-induced poly-
morphic transition should occur, and, if so, at what conditions. The least computationally
demanding is the analysis of a common tangent to the two E(V) curves, as it does not
require the phonon density of states calculations. Other methods include the calculation
and comparison of Gibbs free energy or molecular dynamics simulations.
While, in most of the works, isotropic compression has been used to reproduce the
commonly applied experimental conditions, in the reviewed works, successful applications
of anisotropic compression have also been found.
Finally, it is worth noticing how different levels of theory imply the accuracy of the
results. As in other areas of molecular modeling, in the DFT calculations of organic solids
under high pressure, more computationally demanding calculations and more complex
systems usually yield more reliable results, with the cost being computational time. For
example, while the first-principles methods enable us to calculate lattice energies with
accuracy better than 1 kcal/mol, the requirement is to go beyond a pairwise model of
dispersion correction and include non-additive many-body contributions—for example,
by using the MBD method. However, this approach is significantly more computationally
demanding than the use of the TS or Grimme schemes. Moreover, as described above in
the examples of urea and chlorothiazide, sometimes, to observe the polymorphic phase
transition, molecular dynamics simulations are required, as geometry optimization may not
be sufficient due to the entropy-driven transformations. However, the computational cost
demanding than the use of the TS or Grimme schemes. Moreover, as described above in
the examples of urea and chlorothiazide, sometimes, to observe the polymorphic phase
transition, molecular dynamics simulations are required, as geometry optimization may
Int. J. Mol. Sci. 2023, 24, 14155 not be sufficient due to the entropy-driven transformations. However, the computational 34 of 40
cost of MD simulations is much higher than that of geometry optimization. Another
example is the size of the cell being optimized, as sometimes the modeling of a supercell
composed
of of a fewis(or
MD simulations more)
much unitthan
higher cellsthat
provides more optimization.
of geometry accurate results, especially
Another when
example is
the size
the crystal symmetry
of the cell beingis optimized,
changed during optimization.
as sometimes Unfortunately,
the modeling an increase
of a supercell in the
composed
size
of of the
a few (ormodeled
more) unit system
cellsleads to themore
provides undesired
accurateelongation of calculations.
results, especially when the crystal
The most popular programs that are used to perform such
symmetry is changed during optimization. Unfortunately, an increase in the size calculations are plane-
of the
wave DFT codes, such as CASTEP, VASP, or Quantum
modeled system leads to the undesired elongation of calculations. ESPRESSO (Table 2). Fortunately,
mostTheof them
most are free to
popular use for either
programs non-commercial
that are used to performorsuch academic purposes.
calculations Moreover,
are plane-wave
due to
DFT the variety
codes, such asofCASTEP,
published tutorials
VASP, in the form
or Quantum of both (Table
ESPRESSO videos2).
and pdf documents,
Fortunately, most ofit
is feasible
them totostart
are free use performing such calculations
for either non-commercial even without
or academic the Moreover,
purposes. use of introductory
due to the
courses.of published tutorials in the form of both videos and pdf documents, it is feasible to
variety
Therefore, in
start performing conclusion,
such calculationswe even
would like tothe
without encourage researcherscourses.
use of introductory who are yet not
familiar
Therefore, in conclusion, we would like to encourage researchers who are helpful
with such methods to consider their use, as they can be extremely both in
yet not familiar
planning
with such experiments and during
methods to consider theirthe
use,analysis
as theyand
can interpretation of experimental
be extremely helpful results
both in planning
(Figure 8). and during the analysis and interpretation of experimental results (Figure 8).
experiments

Figure 8.8.AAsummary
Figure summaryof the requirements
of the and possible
requirements outcomes
and possible of the DFT
outcomes calculations
of the of molecular
DFT calculations of
solids under
molecular highunder
solids pressure.
high pressure.

Author Contributions:
Author Contributions: Conceptualization,
Conceptualization, E.N.,
E.N., K.M.
K.M. and
and Ł.S.;
Ł.S.; methodology,
methodology, E.N.,E.N., K.M.
K.M. and
and Ł.S.;
Ł.S.;
software, E.N., K.M. and Ł.S.; validation, E.N., K.M. and Ł.S.; formal analysis, E.N., K.M.
software, E.N., K.M. and Ł.S.; validation, E.N., K.M. and Ł.S.; formal analysis, E.N., K.M. and Ł.S.; and Ł.S.;
investigation, E.N.,
investigation, E.N., K.M.
K.M. and
and Ł.S.;
Ł.S.; resources,
resources, E.N.,
E.N., K.M.
K.M. and
and Ł.S.;
Ł.S.; data
data curation,
curation, E.N.,
E.N., K.M.
K.M. and
and Ł.S.;
Ł.S.;
writing—original draft preparation, E.N., K.M. and Ł.S.; writing—review and editing,
writing—original draft preparation, E.N., K.M. and Ł.S.; writing—review and editing, E.N., K.M. and E.N., K.M.
and Ł.S.;
Ł.S.; visualization,
visualization, E.N.,E.N.,
K.M.K.M. and supervision,
and Ł.S.; Ł.S.; supervision,
E.N., E.N.,
K.M. K.M. and project
and Ł.S.; Ł.S.; project administration,
administration, E.N.,
E.N., K.M. and Ł.S. All authors have read and agreed to the published version
K.M. and Ł.S. All authors have read and agreed to the published version of the manuscript. of the manuscript.
Funding: This research received no external funding.
Funding: This research received no external funding.
Conflicts of Interest: The authors declare no conflicts of interest.
Conflicts of Interest: The authors declare no conflict of interest.

Abbreviations
Abbreviations
Abbreviation Description
2D PES Two-Dimensional Potential Energy Surface
AIM Atoms in Molecules
AS Absorption Spectra Calculations
aiMD Ab Initio Molecular Dynamics
BG Bandgap
Int. J. Mol. Sci. 2023, 24, 14155 35 of 40

CL-20 Hexanitrohexaazaisowurtzitane
CM Center-of-Mass Fractional Position Calculations
COT 1,3,5,7-Cyclooctatetraene
CSP Crystal Structure Prediction
DFT Density Functional Theory
DOS Density of States
EDAB Ethylenediamine Bisborane
EOS Equation of State
ES Excited State Calculation
FOX-7 1,1-Diamino-2,2-Dinitroethene
FPMD First-Principles Molecular Dynamics
GD Grimme Dispersion
GGA Generalized Gradient Approximations
GO Geometry Optimization
HF Hartree–Fock
HMX Cyclic Aliphatic Nitramine Octahydro-1,3,5,7-Tetranitro-1,3,5,7-Tetrazocin
HS Hirshfeld Surface
IGM Intramolecular Gradients Method
INS Inelastic Neutron Scattering
IR Infrared
KS Kohn–Sham
LLM-105 2,6-Diamino-3,5-Dinitropyrazine-1-Oxide
MA Methylammonium
MBD Many-Body Dispersion
MD Molecular dynamics
MO Molecular Orbitals
MPD Mutual Penetration Distance
MSST Multi-Scale Shock Technique
NA Not Applicable
NBO Natural Bond Orbitals
NMR Nuclear Magnetic Resonance
NP Not Provided
NPT Isothermal–Isobaric Ensemble
NVT Canonical Ensemble
OP Optical Properties
PC Phonon DOS Calculation
PD Phase Diagram
PF Phonon Frequency
PL Photoluminescence
pV Pressure–Volume Terms
PXRD Powder X-Ray Diffraction
QHA Quasi Harmonic Approximation
QMD Quantum Molecular Dynamics
RDX 3,5-Trinitrohexahydro-S-Triazine
RMSD Root Mean Square Deviation
SCC—DFTB Self-Consistent Charge Density Functional Tight Binding
SCXRD Single-Crystal X-Ray Diffraction
SOC Spin–Orbit Coupling
SP Single-Point Calculations
TATB Triaminotrinitrobenzene
TB Transition Barrier Calculation
TD Thermodynamics
TD-DFT Time-Dependent Density Functional Theory Calculations
TNBA 2,4,6-Trinitro-3-Bromoanisole
TNT 2,4,6-Trinitrotoluene
TS Tkatchenko–Scheffler
USPEX Universal Structure Predictor: Evolutionary Xtallography
VTST Variational Transition-State Theory
XRD X-Ray Diffraction
ZPVE Zero-Point Vibrational Energy
Int. J. Mol. Sci. 2023, 24, 14155 36 of 40

References
1. Moggach, S.; Parsons, S. Molecular Solids at Extreme Pressure. CrystEngComm 2010, 12, 2515. [CrossRef]
2. Boldyreva, E. High-Pressure Polymorphs of Molecular Solids: When Are They Formed, and When Are They Not? Some Examples
of the Role of Kinetic Control. Cryst. Growth Des. 2007, 7, 1662–1668. [CrossRef]
3. Hoja, J.; Reilly, A.M.; Tkatchenko, A. First-Principles Modeling of Molecular Crystals: Structures and Stabilities, Temperature and
Pressure. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2017, 7, e1294. [CrossRef]
4. Mazurek, A.H.; Szeleszczuk, Ł.; Pisklak, D.M. Periodic DFT Calculations—Review of Applications in the Pharmaceutical Sciences.
Pharmaceutics 2020, 12, 415. [CrossRef]
5. Wu, J.; Yu, D.; Liu, S.; Rong, C.; Zhong, A.; Chattaraj, P.K.; Liu, S. Is It Possible To Determine Oxidation States for Atoms in
Molecules Using Density-Based Quantities? An Information-Theoretic Approach and Conceptual Density Functional Theory
Study. J. Phys. Chem. A 2019, 123, 6751–6760. [CrossRef] [PubMed]
6. Wang, K.; He, X.; Rong, C.; Zhong, A.; Liu, S.; Zhao, D. On the Origin and Nature of Internal Methyl Rotation Barriers: An
Information-Theoretic Approach Study. Theor. Chem. Acc. 2022, 141, 68. [CrossRef]
7. Luo, C.; He, X.; Zhong, A.; Liu, S.; Zhao, D. What Dictates Alkane Isomerization? A Combined Density Functional Theory and
Information-Theoretic Approach Study. Theor. Chem. Acc. 2023, 142, 78. [CrossRef]
8. Hohenberg, P.; Kohn, W. Inhomogeneous Electron Gas. Phys. Rev. 1964, 136, B864–B871. [CrossRef]
9. Kohn, W.; Sham, L.J. Self-Consistent Equations Including Exchange and Correlation Effects. Phys. Rev. 1965, 140, A1133–A1138.
[CrossRef]
10. van Mourik, T.; Bühl, M.; Gaigeot, M.-P. Density Functional Theory across Chemistry, Physics and Biology. Philos. Trans. R. Soc. A
Math. Phys. Eng. Sci. 2014, 372, 20120488. [CrossRef]
11. Perdew, J.P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865–3868.
[CrossRef] [PubMed]
12. Nicholas, A.D.; Zhao, J.; Slebodnick, C.; Ross, N.L.; Cahill, C.L. High-Pressure Structural and Optical Property Evolution of a
Hybrid Indium Halide Perovskite. J. Solid. State Chem. 2021, 300, 122262. [CrossRef]
13. Liu, X.; Michalchuk, A.A.L.; Bhattacharya, B.; Yasuda, N.; Emmerling, F.; Pulham, C.R. High-Pressure Reversibility in a Plastically
Flexible Coordination Polymer Crystal. Nat. Commun. 2021, 12, 3871. [CrossRef]
14. Yesudhas, S.; Burns, R.; Lavina, B.; Tkachev, S.N.; Sun, J.; Ullrich, C.A.; Guha, S. Coupling of Organic Cation and Inorganic Lattice
in Methylammonium Lead Halide Perovskites: Insights into a Pressure-Induced Isostructural Phase Transition. Phys. Rev. Mater.
2020, 4, 105403. [CrossRef]
15. Jiang, H.; Xue, H.; Wang, L.; Tang, F.; Si, F. Effect of Pressure-Induced Structural Phase Transition on Electronic and Optical
Properties of Perovskite CH3NH3PbI3. Mater. Sci. Semicond. Process 2019, 96, 59–65. [CrossRef]
16. Huang, Y.; Wang, L.; Ma, Z.; Wang, F. Pressure-Induced Band Structure Evolution of Halide Perovskites: A First-Principles
Atomic and Electronic Structure Study. J. Phys. Chem. C 2019, 123, 739–745. [CrossRef]
17. Cai, W.; Bi, T.; Zarifi, N.; Terpstra, T.; Zhang, C.; Verdeny, Z.V.; Zurek, E.; Deemyad, S. Effects of Nonhydrostatic Stress on
Structural and Optoelectronic Properties of Methylammonium Lead Bromide Perovskite. J. Phys. Chem. Lett. 2017, 8, 3457–3465.
[CrossRef]
18. Jiang, S.; Fang, Y.; Li, R.; Xiao, H.; Crowley, J.; Wang, C.; White, T.J.; Goddard, W.A.; Wang, Z.; Baikie, T.; et al. Pressure-Dependent
Polymorphism and Band-Gap Tuning of Methylammonium Lead Iodide Perovskite. Angew. Chem.-Int. Ed. 2016, 55, 6540–6544.
[CrossRef]
19. Valiente, R.; García-Lastra, J.M.; García-Fernández, P.; García-Revilla, S.; Wenger, O.S. Red-to-Yellow Pressure-Induced Phase
Transition in Pt(Bpy)Cl2: Spectroscopic Study Supported by DFT Calculations. Eur. J. Inorg. Chem. 2007, 2007, 5735–5742.
[CrossRef]
20. Childs, C.M.; Steele, B.A.; Grivickas, P.; Zhang, D.; Crowhurst, J.; Kuo, I.W.; Bastea, S.; Clarke, S.M. High-Pressure Investigation
of 2,4,6-Trinitro-3-Bromoanisole (TNBA): Structural Determination and Piezochromism. J. Phys. Chem. C 2022, 126, 1176–1187.
[CrossRef]
21. Wu, X.; Li, Y.; Dong, W.; Zhang, J.-G. Pressure-Induced Phase Transition of a Series of Energetic Pentazolate Anion Salts: A DFT
Study. New J. Chem. 2022, 46, 5653–5662. [CrossRef]
22. Gao, C.; Wang, J.; Zhang, Y.; Su, H.; Xu, Z.; Dai, R.; Wang, Z.; Zhang, Z. Pressure-Induced Phase Transition of β-RDX Single
Crystals. J. Phys. Chem. C 2021, 125, 6418–6426. [CrossRef]
23. Guan, R.; Wang, P.; Song, Y.; Staroverov, V.N. Pressure-Induced Polymorphic Transformations of Ethylenediamine Bisborane. J.
Phys. Chem. C 2021, 125, 18614–18622. [CrossRef]
24. Xu, Z.; Chen, Q.; Li, X.; Wang, J.; Wang, X.; Gao, C.; Dai, R.; Wang, Z.; Huang, S.; Liu, Y.; et al. Electronic Structure of LLM-105
Crystal under High Pressure and Low Temperature. J. Phys. Chem. C 2020, 124, 2399–2405. [CrossRef]
25. Su, Y.; Fan, J.; Zheng, Z.; Zhao, J. Uniaxial Compression Behavior and Spectroscopic Properties of Energetic 1,1-Diamino-2,2-
Dinitroethylene (FOX-7) Crystals from Density Functional Theory Calculations. Prog. Prog. Nat. Sci. Mater. Int. 2019, 29, 329–334.
[CrossRef]
26. Konar, S.; Michalchuk, A.A.L.; Sen, N.; Bull, C.L.; Morrison, C.A.; Pulham, C.R. High-Pressure Study of Two Polymorphs of
2,4,6-Trinitrotoluene Using Neutron Powder Diffraction and Density Functional Theory Methods. J. Phys. Chem. C 2019, 123,
26095–26105. [CrossRef]
Int. J. Mol. Sci. 2023, 24, 14155 37 of 40

27. Zhao, X.; Zhu, W. Possible Pre-Phase Transition of the α-HMX Crystal Observed by the Variation of Hydrogen-Bonding Network
under High Pressures. CrystEngComm 2019, 22, 330–348. [CrossRef]
28. Landerville, A.C.; Crowhurst, J.C.; Grant, C.D.; Zaug, J.M.; Oleynik, I. Experimental and Theoretical Investigation of Pressurede-
pendent Raman Spectra of Triaminotrinitrobenzene (TATB) at High Pressures. In Proceedings of the Conference of the American
Physical Society Topical Group on Shock Compression of Condensed Matter, Tampa Bay, FL, USA, 14–19 June 2015; Volume 1793.
29. Dreger, Z.A.; Stash, A.I.; Yu, Z.-G.Z.G.Z.-G.; Chen, Y.-S.Y.S.Y.-S.; Tao, Y.; Gupta, Y.M. High-Pressure Structural Response of an
Insensitive Energetic Crystal: 1, 1-Diamino-2, 2-Dinitroethene (FOX-7). J. Phys. Chem. C 2016, 120, 27600–27607. [CrossRef]
30. Cawkwell, M.J.; Luscher, D.J.; Addessio, F.L.; Ramos, K.J. Equations of State for the α and γ Polymorphs of Cyclotrimethylene
Trinitramine. J. Appl. Phys. 2016, 119, 185106. [CrossRef]
31. Stavrou, E.; Riad Manaa, M.; Zaug, J.M.; Kuo, I.F.W.; Pagoria, P.F.; Kalkan, B.; Crowhurst, J.C.; Armstrong, M.R. The High
Pressure Structure and Equation of State of 2,6-Diamino-3,5-Dinitropyrazine-1-Oxide (LLM-105) up to 20 GPa: X-Ray Diffraction
Measurements and First Principles Molecular Dynamics Simulations. J. Chem. Phys. 2015, 143, 144506. [CrossRef]
32. Manaa, M.R.; Kuo, I.-F.W.; Fried, L.E. First-Principles High-Pressure Unreacted Equation of State and Heat of Formation of
Crystal 2,6-Diamino-3, 5-Dinitropyrazine-1-Oxide (LLM-105). J. Chem. Phys. 2014, 141, 064702. [CrossRef] [PubMed]
33. Ge, N.-N.; Wei, Y.-K.; Zhao, F.; Chen, X.-R.; Ji, G.-F. Pressure-Induced Metallization of Condensed Phase β-HMX under Shock
Loadings via Molecular Dynamics Simulations in Conjunction with Multi-Scale Shock Technique. J. Mol. Model. 2014, 20, 2350.
[CrossRef] [PubMed]
34. Yedukondalu, N.; Vaitheeswaran, G. Polymorphism and Thermodynamic Ground State of Silver Fulminate Studied from van Der
Waals Density Functional Calculations. J. Chem. Phys. 2014, 140, 224705. [CrossRef] [PubMed]
35. Sorescu, D.C.; Rice, B.M. Theoretical Predictions of Energetic Molecular Crystals at Ambient and Hydrostatic Compression
Conditions Using Dispersion Corrections to Conventional Density Functionals (DFT-D). J. Phys. Chem. C 2010, 114, 6734–6748.
[CrossRef]
36. Tkachev, S.N.; Pravica, M.; Kim, E.; Romano, E.; Week, P.F. High-Pressure Studies of 1,3,5,7-Cyclooctatetraene: Experiment and
Theory. J. Phys. Chem. A 2008, 112, 11501–11507. [CrossRef]
37. Xu, X.-J.; Zhu, W.-H.; Xiao, H.-M. DFT Studies on the Four Polymorphs of Crystalline CL-20 and the Influences of Hydrostatic
Pressure on ε-CL-20 Crystal. J. Phys. Chem. B 2007, 111, 2090–2097. [CrossRef]
38. Ali, I.; Tang, J.; Han, Y.; Wei, Z.; Zhang, Y.; Li, J. A Solid-Solid Phase Transformation of Triclabendazole at High Pressures. Crystals
2022, 12, 300. [CrossRef]
39. Szeleszczuk, Ł.; Mazurek, A.H.; Milcarz, K.; Napiórkowska, E.; Pisklak, D.M. Can We Predict the Isosymmetric Phase Transition ?
Application of DFT Calculations to Study the Pressure Induced Transformation of Chlorothiazide. Int. J. Mol. Sci. 2021, 22, 10100.
[CrossRef]
40. Silva, R.; de Castro, A.J.; da Silva Filho, J.G.; de Sousa, F.F.; Paraguassu, W.; Freire, P.T.C.; Façanha Filho, P.F. Pressure-Induced
Phase Transition in Glycinium Maleate Crystal. Spectrochim. Acta A Mol. Biomol. Spectrosc. 2021, 262, 120076. [CrossRef]
41. Cook, C.; McKinley, J.L.; Beran, G.J.O.O. Modeling the α-and β-Resorcinol Phase Boundary via Combination of Density Functional
Theory and Density Functional Tight-Binding. J. Chem. Phys. 2021, 154, 134109. [CrossRef]
42. Xavier, N.F.; Da Silva, A.M.; Bauerfeldt, G.F. What Rules the Relative Stability of α-, β-, and γ-Glycine Polymorphs? Cryst. Growth
Des. 2020, 20, 4695–4706. [CrossRef]
43. Novelli, G.; Maynard-Casely, H.E.; McIntyre, G.J.; Warren, M.R.; Parsons, S. Effect of High Pressure on the Crystal Structures of
Polymorphs of L-Histidine. Cryst. Growth Des. 2020, 20, 7788–7804. [CrossRef]
44. Novelli, G.; Maynard-Casely, H.E.; McIntyre, G.J.; Warren, M.R.; Parsons, S. Correction to “The Effect of High Pressure on the
Crystal Structures of Polymorphs of L-Histidine”. Cryst. Growth Des. 2022, 22, 937–938. [CrossRef]
45. Mazurek, A.; Szeleszczuk, Ł.; Pisklak, D.M. Can We Predict the Pressure Induced Phase Transition of Urea? Application of
Quantum Molecular Dynamics. Molecules 2020, 25, 2021. [CrossRef] [PubMed]
46. Clarke, S.M.; Steele, B.A.; Kroonblawd, M.P.; Zhang, D.; Kuo, I.-F.W.; Stavrou, E. An Isosymmetric High-Pressure Phase Transition
in α-Glycylglycine: A Combined Experimental and Theoretical Study. J. Phys. Chem. B 2020, 124, 1–10. [CrossRef]
47. Giordano, N.; Beavers, C.M.; Kamenev, K.V.; Marshall, W.G.; Moggach, S.A.; Patterson, S.D.; Teat, S.J.; Warren, J.E.; Wood,
P.A.; Parsons, S. High-Pressure Polymorphism in l-Threonine between Ambient Pressure and 22 GPa. CrystEngComm 2019, 21,
4444–4456. [CrossRef]
48. Prazyan, T.L.; Zhuravlev, Y.N.; Golovko, O.V.; Obolonskaya, O.S. DFT-Study of Pressure-Induced Phase Transition in L-Threonine.
J. Mol. Struct. 2019, 1196, 271–279. [CrossRef]
49. Szeleszczuk, Ł.; Pisklak, D.M.; Zielińska-Pisklak, M. Can We Predict the Structure and Stability of Molecular Crystals under
Increased Pressure? First-Principles Study of Glycine Phase Transitions. J. Comput. Chem. 2018, 39, 1300–1306. [CrossRef]
50. Rychkov, D.A.; Stare, J.; Boldyreva, E. V Pressure-Driven Phase Transition Mechanisms Revealed by Quantum Chemistry:
L-Serine Polymorphs. Phys. Chem. Chem. Phys. 2017, 19, 6671–6676. [CrossRef]
51. Fedorov, A.Y.; Rychkov, D.A.; Losev, E.A.; Zakharov, B.A.; Stare, J.; Boldyreva, E.V. Effect of Pressure on Two Polymorphs of
Tolazamide: Why No Interconversion? CrystEngComm 2017, 19, 2243–2252. [CrossRef]
52. Leblanc, L.M.; Otero-de-La-Roza, A.; Johnson, E.R. Evaluation of Shear-Slip Transitions in Crystalline Aspirin by Density-
Functional Theory. Cryst. Growth Des. 2016, 16, 6867–6873. [CrossRef]
Int. J. Mol. Sci. 2023, 24, 14155 38 of 40

53. Adhikari, K.; Flurchick, K.M.; Valenzano, L. A Hybrid Density Functional Study on the Effects of Pressure on Paracetamol and
Aspirin Polymorphs. Comput. Theor. Chem. 2015, 1062, 90–98. [CrossRef]
54. Druzbicki, K.; Mikuli, E.; Pałka, N.; Zalewski, S.; Ossowska-Chruściel, M.D. Polymorphism of Resorcinol Explored by Comple-
mentary Vibrational Spectroscopy (FT-RS, THz-TDS, INS) and First-Principles Solid-State Computations (Plane-Wave DFT). J.
Phys. Chem. B 2015, 119, 1681–1695. [CrossRef] [PubMed]
55. Sabatini, R.; Küçükbenli, E.; Kolb, B.; Thonhauser, T.; De Gironcoli, S. Structural Evolution of Amino Acid Crystals under Stress
from a Non-Empirical Density Functional. J. Phys. Condens. Matter 2012, 24, 424209. [CrossRef] [PubMed]
56. Wood, P.A.; Francis, D.; Marshall, W.G.; Moggach, S.A.; Parsons, S.; Pidcock, E.; Rohl, A.L. A Study of the High-Pressure
Polymorphs of L-Serine Using Ab Initio Structures and PIXEL Calculations. CrystEngComm 2008, 10, 1154–1166. [CrossRef]
57. Howard, C.M.; Wood, I.G.; Knight, K.S.; Fortes, A.D. Ab Initio Simulations of α- and β-Ammonium Carbamate (NH4·NH2CO2),
and the Thermal Expansivity of Deuterated α-Ammonium Carbamate from 4.2 to 180 K by Neutron Powder Diffraction. Acta
Crystallogr. B Struct. Sci. Cryst. Eng. Mater. 2022, 78, 459–475. [CrossRef]
58. Kurzydłowski, D.; Chumak, T.; Rogoża, J. Phase Stability of Chloroform and Dichloromethane at High Pressure. Crystals 2020, 10, 920.
[CrossRef]
59. Batyrev, I.G. The Modeling the Structure of Croconic and Squaric Acids under Pressure. In Proceedings of the Conference of the
American Physical Society Topical Group on Shock Compression of Condensed Matter, Portland, OR, USA, 16–21 June 2019; Lane,
J.M.D., Germann, T.C., A.M.R.W.R.D.D.Z.J., Eds.; American Institute of Physics Inc.: College Park, MD, USA, 2020; Volume 2272.
60. Sahoo, S.; Ravindran, T.R.; Rajaraman, R.; Srihari, V.; Pandey, K.K.; Chandra, S. Pressure-Induced Amorphization of Diiso-
propylammonium Perchlorate Studied by Raman Spectroscopy and X-Ray Diffraction. J. Phys. Chem. A 2020, 124, 1993–2000.
[CrossRef]
61. Zhou, J.; Lian, J.; Ding, Y.; Gao, Y.; Zhu, W. Pressure-Induced Variation of Structural, Elastic, Vibrational, Thermodynamic
Properties and Hardness of C11N4 Polymorphs. Results Phys. 2019, 14, 102453. [CrossRef]
62. Colmenero, F. Mechanical Properties of Anhydrous Oxalic Acid and Oxalic Acid Dihydrate. Phys. Chem. Chem. Phys. 2019, 21,
2673–2690. [CrossRef]
63. Saraiva, G.D.; Maia, J.R.; Lima, J.A., Jr.; Nogueira, C.E.S.; Freire, P.T.C.; de Sousa, F.F.; Teixeira, A.M.R.; Mendes Filho, J. Pressure
Induced Transformations in Sorbic Acid. Spectrochim. Acta A Mol. Biomol. Spectrosc. 2017, 184, 327–334. [CrossRef]
64. Wong, J.W.L.; Mailman, A.; Lekin, K.; Winter, S.M.; Yong, W.; Zhao, J.; Garimella, S.V.; Tse, J.S.; Secco, R.A.; Desgreniers, S.; et al.
Pressure Induced Phase Transitions and Metallization of a Neutral Radical Conductor. J. Am. Chem. Soc. 2014, 136, 1070–1081.
[CrossRef] [PubMed]
65. Lekin, K.; Leitch, A.A.; Tse, J.S.; Bao, X.; Secco, R.A.; Desgreniers, S.; Ohishi, Y.; Oakley, R.T. A Pressure Induced Structural
Dichotomy in Isostructural Bis-1,2,3-Thiaselenazolyl Radical Dimers. Cryst. Growth Des. 2012, 12, 4676–4684. [CrossRef]
66. Schatschneider, B.; Liang, J.-J.J.J. Accurate Modeling of Molecular Crystal through Dispersion-Corrected Density Functional
Theory (DFT-D) Method. In Materials Research Society Symposium Proceedings; MRS Online Proceedings Library (OPL): Warrendale,
PA, USA, 2011; Volume 1301, pp. 125–130.
67. Clark, S.J.; Segall, M.D.; Pickard, C.J.; Hasnip, P.J.; Probert, M.I.J.; Refson, K.; Payne, M.C. First Principles Methods Using CASTEP.
Z. Krist. Cryst. Mater. 2005, 220, 567–570. [CrossRef]
68. CASTEP. Available online: http://www.castep.org/ (accessed on 10 August 2023).
69. Hunter, S.; Sutinen, T.; Parker, S.F.; Morrison, C.A.; Williamson, D.M.; Thompson, S.; Gould, P.J.; Pulham, C.R. Experimental and
DFT-D Studies of the Molecular Organic Energetic Material RDX. J. Phys. Chem. C 2013, 117, 8062–8071. [CrossRef]
70. VASP-Vienna Ab Initio Simulation Package. Available online: https://www.vasp.at/ (accessed on 10 August 2023).
71. Sun, G.; Kürti, J.; Rajczy, P.; Kertesz, M.; Hafner, J.; Kresse, G. Performance of the Vienna Ab Initio Simulation Package (VASP) in
Chemical Applications. J. Mol. Struct. Theochem. 2003, 624, 37–45. [CrossRef]
72. Quantum Espresso. Available online: https://www.quantum-espresso.org/ (accessed on 10 August 2023).
73. Giannozzi, P.; Andreussi, O.; Brumme, T.; Bunau, O.; Buongiorno Nardelli, M.; Calandra, M.; Car, R.; Cavazzoni, C.;
Ceresoli, D.; Cococcioni, M.; et al. Advanced Capabilities for Materials Modelling with Quantum ESPRESSO. J. Phys.
Condens. Matter 2017, 29, 465901. [CrossRef]
74. Gaussian.Com|Expanding the Limits of Computational Chemistry. Available online: https://gaussian.com/ (accessed on
10 August 2023).
75. Frisch, M.J.; Trucks, G.W.; Schlegel, H.B.; Scuseria, G.E.; Robb, M.A.; Cheeseman, J.R.; Scalmani, G.; Barone, V.; Petersson, G.A.;
Nakatsuji, H.; et al. Gaussian 16 Rev. C.01 2016. Available online: https://gaussian.com/citation/ (accessed on 12 September 2023).
76. Kühne, T.D.; Iannuzzi, M.; Del Ben, M.; Rybkin, V.V.; Seewald, P.; Stein, F.; Laino, T.; Khaliullin, R.Z.; Schütt, O.; Schiffmann, F.; et al.
CP2K: An Electronic Structure and Molecular Dynamics Software Package -Quickstep: Efficient and Accurate Electronic Structure
Calculations. J. Chem. Phys. 2020, 152, 194103. [CrossRef]
77. CP2K. Available online: https://www.cp2k.org/ (accessed on 10 August 2023).
78. CRYSTAL. Available online: https://www.crystal.unito.it/index.html (accessed on 10 August 2023).
79. Erba, A.; Desmarais, J.K.; Casassa, S.; Civalleri, B.; Donà, L.; Bush, I.J.; Searle, B.; Maschio, L.; Edith-Daga, L.; Cossard, A.; et al.
CRYSTAL23: A Program for Computational Solid State Physics and Chemistry. J. Chem. Theory Comput. 2022. [CrossRef]
80. DMol3-BIOVIA. Available online: https://www.3ds.com/products-services/biovia/ (accessed on 10 August 2023).
Int. J. Mol. Sci. 2023, 24, 14155 39 of 40

81. Delley, B. DMol3 DFT Studies: From Molecules and Molecular Environments to Surfaces and Solids. Comput. Mater. Sci. 2000, 17,
122–126. [CrossRef]
82. Andreoni, W.; Curioni, A. New Advances in Chemistry and Materials Science with CPMD and Parallel Computing. Parallel
Comput. 2000, 26, 819–842. [CrossRef]
83. CPMD. Available online: http://www.cpmd.org/ (accessed on 10 August 2023).
84. DFTB+. Available online: https://dftbplus.org/ (accessed on 10 August 2023).
85. Hourahine, B.; Aradi, B.; Blum, V.; Bonafé, F.; Buccheri, A.; Camacho, C.; Cevallos, C.; Deshaye, M.Y.; Dumitric, T.; Dominguez,
A.; et al. DFTB+, a Software Package for Efficient Approximate Density Functional Theory Based Atomistic Simulations. J. Chem.
Phys. 2020, 152, 124101. [CrossRef] [PubMed]
86. Soler, J.M.; Artacho, E.; Gale, J.D.; García, A.; Junquera, J.; Ordejón, P.; Sánchez-Portal, D. The SIESTA Method for Ab Initio
Order-N Materials Simulation. J. Phys. Condens. Matter 2002, 14, 2745–2779. [CrossRef]
87. Siesta. Available online: https://siesta-project.org/siesta/ (accessed on 10 August 2023).
88. WIEN2k. Available online: http://www.wien2k.at/index.html (accessed on 10 August 2023).
89. Blaha, P.; Schwarz, K.; Tran, F.; Laskowski, R.; Madsen, G.K.H.; Marks, L.D. WIEN2k: An APW+lo Program for Calculating the
Properties of Solids. J. Chem. Phys. 2020, 152, 074101. [CrossRef]
90. Bučko, T.; Lebègue, S.; Hafner, J.; Ángyán, J.G. Tkatchenko-Scheffler van Der Waals Correction Method with and without
Self-Consistent Screening Applied to Solids. Phys. Rev. B 2013, 87, 064110. [CrossRef]
91. Grimme, S.; Steinmetz, M. Effects of London Dispersion Correction in Density Functional Theory on the Structures of Organic
Molecules in the Gas Phase. Phys. Chem. Chem. Phys. 2013, 15, 16031. [CrossRef]
92. Xu, P.; Alkan, M.; Gordon, M.S. Many-Body Dispersion. Chem. Rev. 2020, 120, 12343–12356. [CrossRef]
93. Winkler, B.; Milman, V. Density Functional Theory Based Calculations for High Pressure Research. Z. Krist. Cryst. Mater. 2014,
229, 112–122. [CrossRef]
94. Baroni, S.; de Gironcoli, S.; Dal Corso, A.; Giannozzi, P. Phonons and Related Crystal Properties from Density-Functional
Perturbation Theory. Rev. Mod. Phys. 2001, 73, 515–562. [CrossRef]
95. Gelpi, J.; Hospital, A.; Goñi, R.; Orozco, M. Molecular Dynamics Simulations: Advances and Applications. Adv. Appl. Bioinform.
Chem. 2015, 37, 37–47. [CrossRef]
96. Starodumov, I.; Ankudinov, V.; Nizovtseva, I. A Review of Continuous Modeling of Periodic Pattern Formation with Modified
Phase-Field Crystal Models. Eur. Phys. J. Spec. Top. 2022, 231, 1135–1145. [CrossRef]
97. Hasnip, P.J.; Refson, K.; Probert, M.I.J.; Yates, J.R.; Clark, S.J.; Pickard, C.J. Density Functional Theory in the Solid State. Philos.
Trans. R. Soc. A: Math. Phys. Eng. Sci. 2014, 372, 20130270. [CrossRef] [PubMed]
98. Hodgkinson, P. NMR Crystallography of Molecular Organics. Prog. Nucl. Magn. Reson. Spectrosc. 2020, 118, 10–53. [CrossRef]
[PubMed]
99. Rychkov, D.A. A Short Review of Current Computational Concepts for High-Pressure Phase Transition Studies in Molecular
Crystals. Crystals 2020, 10, 81. [CrossRef]
100. Fabbiani, F.P.A.; Pulham, C.R. High-Pressure Studies of Pharmaceutical Compounds and Energetic Materials. Chem. Soc. Rev.
2006, 35, 932. [CrossRef]
101. Holy, J.A. Raman Spectroscopy of Aminated and Ultrafine 1,3,5-Triamino-2,4,6-Trinitrobenzene and PBX 9502 as a Function of
Pressing Pressure. J. Phys. Chem. B 2008, 112, 7489–7498. [CrossRef]
102. Oswald, I.D.H.; Lennie, A.R.; Pulham, C.R.; Shankland, K. High-Pressure Structural Studies of the Pharmaceutical, Chlorothiazide.
CrystEngComm 2010, 12, 2533. [CrossRef]
103. Reilly, A.M.; Tkatchenko, A. Role of Dispersion Interactions in the Polymorphism and Entropic Stabilization of the Aspirin
Crystal. Phys. Rev. Lett. 2014, 113, 055701. [CrossRef]
104. Crowell, E.L.; Dreger, Z.A.; Gupta, Y.M. High-Pressure Polymorphism of Acetylsalicylic Acid (Aspirin): Raman Spectroscopy. J.
Mol. Struct. 2015, 1082, 29–37. [CrossRef]
105. Marom, N.; DiStasio, R.A.; Atalla, V.; Levchenko, S.; Reilly, A.M.; Chelikowsky, J.R.; Leiserowitz, L.; Tkatchenko, A. Many-Body
Dispersion Interactions in Molecular Crystal Polymorphism. Angew. Chem. Int. Ed. 2013, 52, 6629–6632. [CrossRef]
106. Rodríguez, J.S.; Costa, G.; da Silva, M.B.; Silva, B.P.; Honório, L.J.; de Lima-Neto, P.; Santos, R.C.R.; Caetano, E.W.S.; Alves, H.W.L.;
Freire, V.N. Structural and Optoelectronic Properties of the α-, β-, and γ-Glycine Polymorphs and the Glycine Dihydrate Crystal:
A DFT Study. Cryst. Growth Des. 2019, 19, 5204–5217. [CrossRef]
107. Dudek, M.K.; Drużbicki, K. Along the Road to Crystal Structure Prediction (CSP) of Pharmaceutical-like Molecules. CrystEngComm
2022, 24, 1665–1678. [CrossRef]
108. Nyman, J.; Reutzel-Edens, S.M. Crystal Structure Prediction Is Changing from Basic Science to Applied Technology. Faraday
Discuss. 2018, 211, 459–476. [CrossRef] [PubMed]
109. Reilly, A.M.; Cooper, R.I.; Adjiman, C.S.; Bhattacharya, S.; Boese, A.D.; Brandenburg, J.G.; Bygrave, P.J.; Bylsma, R.; Campbell, J.E.;
Car, R.; et al. Report on the Sixth Blind Test of Organic Crystal Structure Prediction Methods. Acta Crystallogr. B Struct. Sci. Cryst.
Eng. Mater. 2016, 72, 439–459. [CrossRef] [PubMed]
110. Conclusions of the 7th Crystal Structure Prediction Blind Test|CCDC. Available online: https://www.ccdc.cam.ac.uk/discover/
news/conclusions-of-the-7th-crystal-structure-prediction-blind-test/ (accessed on 4 September 2023).
Int. J. Mol. Sci. 2023, 24, 14155 40 of 40

111. Bhardwaj, R.M.; McMahon, J.A.; Nyman, J.; Price, L.S.; Konar, S.; Oswald, I.D.H.; Pulham, C.R.; Price, S.L.; Reutzel-Edens, S.M. A
Prolific Solvate Former, Galunisertib, under the Pressure of Crystal Structure Prediction, Produces Ten Diverse Polymorphs. J.
Am. Chem. Soc. 2019, 141, 13887–13897. [CrossRef]
112. Neumann, M.A.; van de Streek, J.; Fabbiani, F.P.A.; Hidber, P.; Grassmann, O. Combined Crystal Structure Prediction and
High-Pressure Crystallization in Rational Pharmaceutical Polymorph Screening. Nat. Commun. 2015, 6, 7793. [CrossRef]
113. Jaffe, A.; Lin, Y.; Mao, W.L.; Karunadasa, H.I. Pressure-Induced Conductivity and Yellow-to-Black Piezochromism in a Layered
Cu–Cl Hybrid Perovskite. J. Am. Chem. Soc. 2015, 137, 1673–1678. [CrossRef]

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual
author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to
people or property resulting from any ideas, methods, instructions or products referred to in the content.

You might also like