Theoretical study on the mechanisms of formation of primal carbon clusters and nanoparticles in space

Abstract

We present a theoretical study of assembling clusters and nanoparticles in space from primordial aggregations of unbound carbon atoms. Geometry optimization and SCC-DFTB dynamics methods are employed to predict carbon clusters, their time evolution and stability. The initial density of the aggregates is found to be of primary importance for the structure of the clusters. Aggregates with low initial density yield clusters with an approximately equal prevalence of sp and sp² hybridization with almost missing sp³. Higher initial density results in sp²-dominant molecules, resembling the carbon skeleton of polycyclic aromatic hydrocarbons (PAHs). Larger initial aggregations result in sp²-dominant polymers. Such materials are highly porous and possess a similarity to laterally bound nanotubes. Some clusters resemble fullerene building blocks. We employed metadynamics to model the inter-fragment coupling of such structures and predict the formation of spheroid nanoparticles, closely resembling fullerenes. One such structure has the lowest binding energy per atom among the studied molecules. All zero-dimensional forms, obtained by the simulations, conform to the experimentally detected types of molecules in space. The theoretical IR spectrum of the nanoparticles closely resembles that of fullerene C₇₀ and therefore such imperfect structures may be mistaken for known fullerenes in experimental infrared (IR) telescope studies.

---

1. Introduction

More than any other element, carbon has the innate ability to catenate into diverse, structurally rich chains of arbitrary length. Unsurprisingly, the possible organic molecules outnumber the inorganic ones by two orders of magnitude. Property-wise, carbon-based chemistry has no rival. In fact, for over a decade, it has been expanding in experimental and industrial fields, typically designated for inorganic materials.¹ Form variety and reactivity under mild conditions enable carbon as the fundamental building block of biological systems – the most complex structures known. A basic fundamental question is that of the formation and reformation of carbon phases in the vast and harsh conditions of space. Astrophysics has been answering questions about the origin and distribution of this intriguing element in the Universe. However, it is up to chemistry to answer what would be the very first structures carbon atoms form upon their first meeting in conditions enabling reactivity towards complex molecules. This is the scope of the current scientific study. The origin of carbon in the Universe lies in the triple-alpha fusion process in stellar cores. Star collapse and supernova explosions disperse the atoms in the interstellar and intergalactic space. While circumstellar envelopes are considered the largest reactors for carbon chemistry in the cosmos,^2,3 pathways to complex carbon molecules also exist in interstellar clouds.^4,5 Carbonaceous compounds are also found in distant galaxies^6,7 and protoplanetary nebulae.^8–10 Although the detection of carbon in exoplanetary atmospheres has many obstacles, and some claims are being debated,¹¹ some are considered reliable.^12–14 Surprisingly, polycyclic aromatic hydrocarbons (PAHs) were found even in the cold, dark Taurus Molecular Cloud 1 (TMC-1), at a temperature of 10 K, far away from any evolved star.¹⁵ In our solar system, carbon has been registered in the atmospheres of planets and moons,^16–18 and even in comets.^19–21

To date, of 256 molecules detected in the cosmos, the existence of 204 is proven, and most of them contain carbon.^22,23 The size of the structures ranges from 2 to 70 atoms.²³ The first carbon chemical to be discovered in space is CH (methylidyne);²⁴ CN (cyan radical),²⁵ CS (carbon monosulfide)²⁶ and C₂²⁷ are other examples of such small molecules. Clusters of large carbon content can give us an insight into the preferred structures in cosmic environments consisting primarily of the element of life: L-C₃H⁺ (cyclopropynylidynium cation),²⁸ HCCO (ketenyl radical),²⁹ C₄H (butadiynyl radical),³⁰ C₅ (pentacarbon),³¹ H₂CCC (propadienylidene),³² C₅N (cyanobutadiynyl radical),³³ C₆H (hexatriynyl radical),³⁴ HC₇N (cyanotriacetylene),³⁵ C₈H (octatriynyl radical),³⁶ HC₉N (cyanotetraacetylene),³⁷ C₆H₆ (benzene)³⁸ and C₆H₅CN (benzonitrile).³⁹ Intriguing small to medium-sized cyclic molecules have been detected in TMC-1: C₃H₂ (cyclopropenylidene),⁴⁰ C₅H₆ (cyclopentadiene),⁴⁰ C₅H₅CN (cyanocyclopentadiene),⁴¹ C₅H₅CCH (ethynyl cyclopentadiene)⁴² and o-C₆H₄ (ortho-benzyne).⁴³ Even more intriguing are the polycyclic species detected in this nebula: C₉H₈ (indene)⁴⁰ and two isomers of C₁₁H₇N (1- and 2-cyanonaphthalene).¹⁵ All experimentally known allotropes of carbon (diamond, graphite, and fullerenes) have been detected in space.⁴⁴ Of special interest for structural variety and complexity is the discovery of fullerenes in the cosmos. C₆₀ has been detected in the interstellar space^45,46 and in planetary and protoplanetary nebulae.^45,47,48 Its cation also exists in the diffuse interstellar medium.⁴⁵ C₇₀ has been found in a planetary nebula in Tc1.⁴⁵ Multiple sources of rugbyballene's IR signal have also been detected in planetary nebulae of our satellite galaxies – The Magellanic Clouds.⁴⁹

Since experimentally simulating the exact harsh environment of cosmic space and specifics in the events of original carbon atom aggregation poses many obstacles, we turn to theoretical chemistry to answer the postulated question of initial carbon clusters. Quantum chemistry has a decades-proven track record of modeling reactivity and explaining the relationships between structure and properties. Carbon-wise computerized research has been successfully employed in applications ranging anything from predicting and modulating physical and chemical properties of novel theoretical allotropes,^50–52 through phase transitions of solid-state materials,⁵³ to complete ground state and excited state reactivity of organic molecules, and reactions mechanism characterization.⁵⁴

Models of carbon allotropes, such as D-carbon,⁵⁵ modulated T-carbon,⁵⁶ and ultra-hard rhombohedral carbon,⁵⁷ have been employed to yield total energies, pressure stability, and band gaps. Datta et al.⁵⁸ have investigated the effects of the introduction of Stone–Wales defects on the metallicity of carbon nanotubes. Marchant et al.⁵⁹ have explored the configurational space of elemental carbon to determine macroscopic properties, melting transitions, and thermodynamically stable structures across a wide range of pressures. A simple simulation code for the formation of planar hydrocarbon clusters has been proposed.⁶⁰ The IR spectra of a generated population of C₆₀ isomers have been compared to that of buckminsterfullerene, to discover that the plateau in the 6–9 μm region can only be observed from closed-cage structures.⁶¹ The optical spectra of generations of C_n=24,42,60 clusters have been assigned over structural characteristics, such as asphericity and hybridization, to provide an alternative explanation for a UV bump in the interstellar medium extinction curve of galaxies.⁶² The possibility of formation of amorphous hydrocarbon nanoparticles from PAHs has been proposed.⁶³ Amorphous carbon at low densities has been evaluated in terms of hybridization prevalence, density of vibrational states, and localized vibrational modes⁶⁴ – this appears to be the only previous research using disordered atom arrangements as starting structures, instead of accepted experimental or theoretical carbon formations. Another theoretical study investigates the relationship between density and structural characteristics of amorphous carbon, such as hybridization prevalence and interatomic orientation.⁶⁵ High-level multi-configurational theory has revealed that the cyclic isomer of C₃⁺ is 5.2 kcal mol⁻¹ more stable than the linear one.⁶⁶ Takahashi⁶⁷ has predicted certain small-sized C_n (n = 2–8) interstellar carbon species, focusing on stability and proposing mechanisms of formation. The 3-atomic ring systems have been found more stable than the linear isomers. A theoretical model has been proposed for the mechanism of formation of carbon nanotubes in the cosmos.⁶⁸ Two major approaches exist, as propositions on the formation of fullerenes in space. The “top–down” approach starts with the highly ordered graphene,^69–71 while the “bottom-up” approach initiates with smaller molecules.^72–75 A commentary article from 2022 stresses the open question of formation and destruction mechanisms for carbon molecules in the cosmos, as well as challenges in spectral band assignment.⁷⁶ To the present day, the theoretical investigation of carbon clusterization is incomplete, hence the importance of the current research lies in the fulfilment of the task.

This study is focused on the theoretical investigation and characterization of the paths of formation of original carbon clusters in space. The characterization of the resulting systems involves the role of orientation, size, and density of the initial, pre-reaction atomic aggregations. The cluster description includes the prevalence of structural elements and hybridization. The stabilization of the systems was conducted using DFTB molecular dynamics. The resulting molecules, resembling building blocks for closed-cage structures, were assembled into semi-spherical nanoparticles with DFTB metadynamics. To our knowledge, this state-of-the-art method for reaction mechanism research has previously not been employed for carbon clusterization. The theoretical IR analysis indicates that even semi-spherical nanoparticles, assembled from disordered initial atomic formations, can be misidentified by experimental means as known fullerenes. Although the primary focus is carbon chemistry in space, the modeled reactions can also occur in other mediums, such as combustion, plasma cooling and similar high-energy or low pressure environments. According to our literature survey, this is the first theoretical model proposing and investigating the assembly of fullerene-like structures from primordial, unbound aggregations of carbon atoms.

2. Methodology

All calculations are performed using the CP2K/Quickstep package.^77,78 The SCF optimizations are completed using the self-consistent charge density functional based tight binding (SCC-DFTB/DFTB2) method.⁷⁹ An efficient a posteriori treatment for dispersion interactions is employed.⁸⁰

The DFTB method is an approximation to density functional theory (DFT), where the Kohn–Sham (KS) equations are transformed into a form of tight binding equations⁸¹ related to the Harris functional.⁸² The second-order expansion of the KS equations used here enables a transparent, parameter-free generalized Hamiltonian matrix. Its elements are modified by a self-consistent redistribution of Mulliken charges (SCC).⁸³ The KS energy additionally includes a Coulomb interaction between charge fluctuations. The accuracy of DFTB with SCC expansion is comparable to that of DFT methods and higher levels of ab initio theory for various properties of single molecules, solutions, and solid-state materials. The method yields satisfactory geometries and total energies.^84–86 DFTB produces good charge distribution, binding energies, and vibrational frequencies of charged solvated species.⁸⁷ Activation energies in organic chemistry also conform to higher levels of theory,⁸² enabling the study of reaction mechanisms. It has been found that the method yields excellent geometries and energetics for pure carbon species, such as fullerenes ranging from C₂₀ to C₈₆.⁸⁸ At the same time, the method is orders of magnitude less computationally intensive than DFT. The binding energies per C atom for all C_n clusters formed are calculated as E_b = [E(C_n) − n × E(C)]/n, where E(C_n) and E(C) are the energies of the corresponding cluster and an isolated carbon atom, respectively, while n is the number of C atoms in the cluster. With this definition, a lower (more negative) E_b value corresponds to a more stable C_n cluster.

All simulations are performed using Born–Oppenheimer molecular dynamics (BOMD).⁸⁹ The method of metadynamics (MTD) is chosen to model the formation of nanoparticles.^90,91 Metadynamics is a state-of-the-art simulation method in which the processes are guided to model a chosen chemical reaction. Collective variables (colvars, CVs) are defined over molecular degrees of freedom to bias the system towards selected changes. Penalty potentials (hills) are periodically spawned for current values in the CV space to raise the free energy and reach unexamined geometries. When a TS is crossed over, the study of the new minimum begins, once again starting at the bottom of the energy well. Reversing the bias potential peaks gives us the relative stability of each geometry, effectively mapping the free energy surface of the reaction. With the rise in energy, unguided changes are free to occur, giving a realistic insight into the studied processes.

Metadynamics simulations are carried out in the NVT ensemble, using a canonical sampling through velocity rescaling (CSVR) thermostat,^92,93 set to a temperature of 400 K. The timestep is 1 fs. The height of the Gaussian penalty potentials is 1.255 kcal mol⁻¹. The scale factor (Gaussian width) for each collective variable is 0.2. Hills are spawned every 50 fs. All walls are of quadratic type with a potential constant of 20 kcal mol⁻¹. The temperature tolerance is always set to 50 K. The shortest intermolecular distances are above 3.4 Å in all initial (zeroth) steps. Each MTD run is preceded by the geometry optimization of the systems.

All modelled molecules are in singlet spin state.

All IR spectrum results are from standard, static Hessian calculations.

3. Results and discussion

3.1 Structural tendencies for carbon clusters in space

The structural prevalence of carbon clusters is studied by a theoretical examination of different initial atom aggregations. Each such formation undergoes two computational steps to model the natural appearance of primal carbon molecules. The first step is geometry optimization. The second step is BOMD, enabling spontaneous reformation into a more stable form. Dynamical simulations are carried out in the NVT ensemble if the goal is a lone cluster, or in the NPT ensemble to enable polymer formation. The NPT ensemble tailors the cell dimensions to those of the modeled system and allows intercellular covalent bonding within multidimensional systems. Atom hybridization is determined on the basis of the number of neighboring atoms, orbital orientation, bond lengths and three-point angles. The binding energies per C atom (E_b) are estimated after geometry optimization of the final dynamics structures.

3.1.1 Carbon cluster I. The first studied system consists of 25 carbon atoms, randomly positioned in 3D space and centered in a 15.0 × 15.0 × 15.0 Å cell (Fig. 1a). An interpreter language script with a random number generator is used to yield the Cartesian coordinates of each atom in the initial geometry of the aggregate. The original density is 1.25 g cm⁻³. The model is of isolated atoms in vacuo. Geometry optimization yields a cluster with approximately 2D geometry, resembling the C-skeleton of a polycyclic aromatic hydrocarbon (PAH). According to molecular topology and bond lengths, sp² hybridization is around three times more dominant than sp³. There are no atoms in the sp state. The molecule consists of networks of conjugated double bonds. To ensure the stability of the final structure (Fig. 1c), 20 ps of DFTB molecular dynamics in the NVT ensemble is carried out. A single reaction with a minor structural significance occurs: the breaking of a C–C bond. The product, named CCI, conforms to types of cyclic molecules, experimentally detected in space,^{38,39,41–43} especially the polycyclic indene.⁴⁰ The E_b is estimated to be −7.52 eV. Additionally, CCI has a resemblance to a fullerene building block and is used for the assembly of a closed-cage structure (later in the article).

Fig. 1 Geometries of the initial carbon atom arrangement (a) and (d); structures after geometry optimization (b) and (e); and the final forms after DFTB molecular dynamics (c) and (f) for CCI (a)–(c), and CCIITMP (d)–(f). The two isomers at the bottom are to visualize the reformation of MTD of CCIITMP (g) to CCII (h): the distance CVs are shown in violet and the unbiased chemical changes are shown in blue.

3.1.2 Carbon cluster II. The second studied system consists of 22 carbon atoms, arranged in a tetrahedral pattern and centered in a 15.0 × 15.0 × 15.0 Å cell (Fig. 1d). The initial density is 1.37 g cm⁻³. The model is of isolated atoms in vacuo. Geometry optimization results in a tangled polycyclic 3D network of rings with varying sizes (Fig. 1e). Molecular topology and bond lengths indicate that atom hybridization has the following order of prevalence: sp² > sp³ > sp. The carbon atoms in the sp² state are almost twice as those in the sp³ state. During 20 ps of DFTB molecular dynamics, in the NVT ensemble, ring-opening and ring-closing reactions yield significant structural rearrangements. The resultant molecule (CCIITMP) resembles the C-skeleton of a PAH and its geometry (Fig. 1f) is close to planar. The predominant hybridization of the carbon atoms remains sp². CCIITMP somewhat resembles a fullerene fragment. Metadynamics is employed to reform the structure into a building block for a closed-cage molecule, closer to experimentally known fullerenes. Two distance colvars (Fig. 1g) are set. In the course of the simulation, unbiased chemical processes (Fig. 1g) result in a cluster of satisfactory geometry (named CCII, Fig. 1h). Only one of the anticipating covalent changes according to the CVs was necessary. The free energy profile of the occurred guided reaction is shown in Fig. 6a. The process is exothermic, with a product quite more stable than the reagent: the forward barrier is 9 kcal mol⁻¹, while the reverse barrier is 32 kcal mol⁻¹. The E_b of CCII is estimated to be −7.63 eV. This molecule is a single conjugated system of double and triple bonds. The cluster is used for the assembly of a closed-cage, fullerene-like structure (later in the article).

3.1.3 Carbon cluster III. The initial formation consists of 20 carbon atoms, arranged into a somewhat symmetrical geometry (Fig. 2a). The system is centered in a 16.0 × 16.0 × 16.0 Å cell. The mean interatomic distance is approximately twice that of the aggregation to result in CCI. The initial density is 0.57 g cm⁻³. Geometry optimization results in the structure shown in Fig. 2b. According to topology and bond lengths, half of the atoms are in sp hybridization and the other half are in sp². There are no carbons in the sp³ state. Two rings are formed and both are 3-atomic. The structure is mostly branched linear. DFTB molecular dynamics in the NVT ensemble yields no chemical reactions in 20 ps (Fig. 2b). E_b is found to be −7.15 eV. The linear fragments of the structure correspond to the experimentally detected molecules in space.^31–37 According to our results, it is expected that such molecules appear when the original density of the carbon aggregation is low. The spontaneous formation of 3-atomic rings conforms to a previous theoretical conclusion that such rings are energetically preferred to linear structures when it comes to chains of up to 7 carbon atoms.⁶⁷

Fig. 2 Geometries of the initial carbon atom arrangement (a) and (c); structures after geometry optimization (d) and the final cluster (b) and (e) after DFTB molecular dynamics for CCIII (a) and (b), and CCIV (c)–(e). Two 3 × 1 × 3 multicell views of the 2D polymer CCIV: (f) and (g).

3.1.4 Carbon cluster IV. A larger cluster of approximately double the previous sizes is simulated in a supercell with 3D periodic boundary conditions (PBCs). The initial geometry is 40 carbon atoms in a tetrahedral arrangement (Fig. 2c), centered in a 13.5 × 16.5 × 13.5 Å cell. The model is of isolated atoms in vacuo. The initial density is higher than that of CCIII and close to the densities of CCI and CCII (1.6 g cm⁻³). Geometry optimization yields the structure shown in Fig. 2d. The structure is tangled polycyclic. The order of prevalence of hybridization is sp³ > sp² > sp. There are almost as many atoms in the sp² state as in sp³. In our entire research, this is the only structure in which the dominant hybridization is sp³.

DFTB molecular dynamics is carried out for 20 ps in an NPT ensemble. A supercell with initial dimensions of 9.0 × 12.0 × 10.0 Å is used with 3D PBCs. After many ring-opening and ring-closing reactions, the structure is significantly reformed into a 2D polymer. The elementary cell is a complex network of fused rings, arranged in a somewhat tubular formation (Fig. 2e). The single cell dimensions are 8.02 × 10.69 × 8.91 Å. The predominant hybridization is once again sp². E_b is estimated to be −8.17 eV. The polymer somewhat resembles a network of laterally bound carbon nanotubes with a weak-defined structure (Fig. 2f and g).

The event of carbon atoms congregating in vacuo results in various formations, depending on the number of atoms and the initial density. Lone carbon aggregations of 20 to 25 atoms with sufficiently high-density form sp²-dominant clusters representing the C-skeleton of PAHs. This result is independent of the initial interatomic orientation (for example, random or tetrahedral). Such clusters have resemblance to fullerene building blocks. Similar-sized aggregations at sufficiently reduced density result in branched molecules with far more prominent sp hybridization, a very similar degree of sp² hybridization, no sp³ hybridization, and very few rings of low atomic number. If the PBC model has varying cell dimensions, tailored by the cluster size, sp²-dominant multi-dimensional polymers are formed at the initial higher density. Such materials exhibit a porous structure and somewhat resemble laterally bound nanotubes. It is possible that at really low temperatures and/or short periods, sp³-dominant tangled-polycyclic clusters occur in cases of large, pre-reaction carbon aggregations, as evident from the pre-dynamics optimization of CCIV. DFTB molecular dynamics runs sometimes lead to the stabilization and reformation of the structures. The procedure tends to increase the number of atoms in sp² hybridization and decrease the number of atoms in sp³ hybridization, if originally present. Hence, we also expect such tendency to guide spontaneous transformations in nature. Structures arising from initial aggregations with higher density tend to be more stable, as indicated by the binding energy per C atom. Enabling polymerization, with a theoretically unlimited number of atoms, amplifies this tendency. Spontaneous covalent coupling between smaller clusters is expected, in the event of their gathering. All 0-dimensional clusters conform to experimentally detected types of carbon molecules in space.

3.2 The formation of spherical carbon nanoparticles

The binding of clusters into nanoparticles is simulated using metadynamics.

3.2.1 Modeling the synthesis of nanoparticle C₇₁SC. Metadynamics is employed in an attempt to react identical CCI instances into a spherical nanoparticle. At first 2 CCI clusters (Fig. 3a) are bound to a half-sphere. The fragment-coupling reactions are guided with four C–C distance CVs. Three of the interfragment bonds appear almost without a barrier (<2 kcal mol⁻¹, Fig. 4a–c), the fourth interfragment bond requires 11 kcal mol⁻¹ due to the need for a small whole-system twist (Fig. 4d). In the reactions guided with CV2, a neighboring atom participates in the formation of a bond instead of a targeted atom, and an artifact appears in the free energy profile – the covalent interatomic distance minimum has an unexpectedly large value of 2.5 Å. Since penalty potentials are automatically generated at even intervals, perhaps the negligible barriers found are simply an artifact of the time required for the clusters to bind at product bond lengths. A product is formed in less than 1 ps. The resulting cluster, named simply 2 × CCI (Fig. 3b), has 34 atoms in its round fragment and 16 atoms in side chains. The product free energy well is not completely studied, but a reverse barrier of 37 kcal mol⁻¹ at the end of the simulation guarantees that the product is an order of magnitude energetically more stable than the reagents.

Fig. 3 Major structures in the mechanism of formation of the semi-spherical, closed cage carbon nanoparticle C₇₂SC: two CCI fragments (a), 2 × CCI (b), two 2 × CCI fragments (c), C₇₁SC (d) and C₇₂SC (e).

Fig. 4 Free energy profiles of the interfragment coupling reactions in the formation of 2 × CCI: (a)–(d); the formation of C₇₁SC: (e)–(h); and its isomerization to C₇₂SC (i)–(l). All CV values are in units of Å.

Two 2 × CCI units (Fig. 3c) react with each other in another MTD step to a nanoparticle (NP) with an entirely spherical fragment. The fragment-coupling reactions are guided with four C–C distance CVs. Three of the interfragment bonds appear in an almost barrierless fashion and a product is formed in less than 1 ps. The small barriers found (<2 kcal mol⁻¹, Fig. 4e, Fig. 4g and h) are once again probably an artifact of the time required for the fragments to bind at product bond lengths. After the first three of the bonds are formed, the positions of the two carbon atoms in CV2 are less appropriate for coupling which leads to a single barrier of 11 kcal mol⁻¹ (Fig. 4f). Additional fragment–binding bonds are formed without CV bias, resulting in a structure with 71 atoms with a spherical shape (C₇₁SC, Fig. 3d). The barrier for the reverse process is found to be 37 kcal mol⁻¹, which shows the notable stability of the C₇₁SC structure.

3.2.2 Structural improvements of C₇₁SC. Improvements in the structure of C₇₁SC are carried out using MTD to achieve a more symmetrical shape. A few gaps are to be closed by intramolecular bond formation. Four C–C distance CVs are used and the barriers found are: 0.5, 2.5, 7, and 24 kcal mol⁻¹ (Fig. 4k, l, j and i). The product-free energy well is once again not entirely explored but the barrier of 32 kcal mol⁻¹ at the end of the simulation for the reverse process means that the product is significantly more stable than the reagent and gap-reopening reactions are very unlikely. The final fullerene-like nanoparticle is C₇₂SC and represents a C₇₂[4,5,6,7,8,10]fullerene with side chains (Fig. 3e). E_b is found to be −8.01 eV.

We consider our model for the synthesis of closed-shell particles in outer space as a representation of possible natural processes for the following reasons: (i) while almost all the interstellar medium can be characterized with extremely low temperatures and densities, we aim to simulate environments in which rich chemistry may occur: furthermost areas of star coronas, planets, planetary nebulae and warm (lit/irradiated) molecular clouds; (ii) geometrical optimizations yield a fairly probable result for the “collapse” of non-covalent aggregates of carbon atoms, as the major guiding factor in such computation is the energy descent; (iii) BOMD simulations allow no artificial interference in the natural evolution of the chemical systems and entirely follow the Newtonian laws; (iv) the reaction bias, which transforms pure dynamics into MTD, serves only to accelerate a selected reaction and does not exclude alternative processes. Although the rise of energy, due to the creation of penalty potentials (distributed by the CVs), is aimed at chosen chemical events, if a parallel/conjunctive process requires lower energy to occur – it will occur during the simulation. The total energy of the molecular system is increased, solely to facilitate the crossing of a transition state and to explore a new energy well using the same model. MTD is simply an accelerated version of standard dynamics, so the chances for chemical reactions do not remain purely statistical and the simulation may finish in a feasible amount of time; (v) the CVs, chosen for assembly between primal carbon clusters, are only distances between outer-most, least saturated (in terms of total valence) carbon atoms. The energy barriers of the reactions with such atoms are the lowest, as C≡C species are the most reactive. Partial radical character is also expected for such carbon atoms, which would reduce the barriers further. The reactivity of edge atoms causes the least amount of changes in the already established structures, once again pertaining to energy barriers and the probability of process occurrence. Furthermore, the selected atoms are open for more angles of attack; (vi) apt to the aforementioned expectations, all modeled reactions in the assembly of closed-cage structures are barrierless; (vii) the existence of complex, pure carbon molecules in the cosmos proves that carbon-only non-covalent aggregates do occur in space; (viii) as modeled in this article, low density, carbon aggregates do not stabilize into molecules, which resemble building blocks for closed-cage structures. Additionally, the atom-by-atom growth of the existing carbon clusters would meet at least two challenges. The first is that single carbon atoms are too reactive to attach only to edge atoms. The second is that the most stable final configuration for edge-addition is three-membered rings; (xix) points (vii) and (viii) lead to the conclusion that events of spontaneous aggregation of carbon atoms and the following processes of stabilization into small clusters (such as C₂₅ and C₂₂) are a probable pathway to the assembly of closed-cage structures; (x) entropy factors in space enable the occasional proper orientation between primal carbon clusters, enabling the assembly of fullerene-like nanoparticles.

3.2.3 Modeling the synthesis of nanoparticle C₆₆. Another spherical nanoparticle is formed from CCII clusters with the use of metadynamics. The first step is the binding of two CCII fragments (Fig. 5a). The fragment-coupling reactions are guided with three C–C distance CVs. All three reactions are almost barrierless (Fig. 6b–d). The small barriers found (<2 kcal mol⁻¹) are once again probably an artifact of the time required for the fragments to bond at product bond lengths. A fourth C–C coupling reaction occurs without CV bias. All three set reactions are exothermic. The resulting cluster is formed in less than 1 ps of simulation. Although the product free energy well is not completely explored, according to the final energy profiles, the product is an order of magnitude more stable than the reagents. The product is a half-sphere named 2 × CCII (Fig. 5b). This molecule is a single chain of 44 atoms arranged in a bowl formation.

Fig. 5 Major structures in the mechanism of formation of the semi-spherical, closed cage carbon nanoparticle C₆₆: two CCII fragments (a), 2 × CCII (b), one CCII with one 2 × CCI fragment (c), 3 × CCII (d) and C₆₆ (e).

Fig. 6 Free energy profile of the MTD isomerization of CCIITMP to CCII (a). Free energy profiles of the interfragment coupling reactions in the formation of 2 × CCII: (b)–(d); and the formation of 3 × CCII: (e)–(h). All CV values are in units of Å.

An additional CCII reacts with 2 × CCII (Fig. 5c) in another MTD simulation to form a semi-spherical NP. The fragment-coupling reactions are guided with four C–C distance CVs. All four reactions are almost barrierless (Fig. 6e–h). A product is formed in less than 1 ps. Further fragment-coupling bonds are formed without CV bias, resulting in a spherical structure (3 × CCII, Fig. 5d) with 66 atoms and no side chains. The barrier for the reverse process is as high as 27 kcal mol⁻¹, manifesting the significant stability of the formed 3 × CCII structure.

3.2.4 Structural improvements of 3 × CCII. Improvement in the structure of 3 × CCII is carried out using MTD for a more symmetrical shape. A few gaps are to be closed by intramolecular bond formation. A total of three C–C distance CVs are used. The largest forward free energy barrier found is 4 kcal mol⁻¹ (Fig. 7) and gap-reopening reactions are very unlikely. The final fullerene-like nanoparticle is C₆₆ and represents a C₆₆[4,5,6,7,8,13]fullerene (Fig. 5e). Its binding energy per atom is the lowest among the studied molecules at −8.26 eV. The cluster is not as symmetrical as fullerenes C₆₀ or C₇₀ but considering the fact that the initial state of aggregation involved carbon atoms scattered in a tetrahedral orientation to each other and the final result is significant.

Fig. 7 Free energy profiles of the gap-closing reactions for the isomerization of 3 × CCII to C₆₆. All CV values are in units of Å.

In various “top–down” studies on fullerene synthesis, the pre-final step involves irregular, closed-cage structures like C₇₂SC and C₆₆.^66–68 With our study, we correct a previous hypothesis over the mechanism of assembly of fullerenes in space:⁹⁴ (1) graphene is not a required intermediate and (2) tangled and flat polycyclic structures do not reform only in a direction away from fullerenes. A currently accepted idea is that molecules, similar to the well-known fullerenes, eventually transform in space into C₆₀ and C₇₀, etc. These transformations occur because the later can be considered as islands of stability among the closed-cage carbon structures.⁹⁴

3.3 IR-spectrum resemblance to known fullerenes

IR spectroscopy is the primary means of detecting fullerenes in space. Multiple sources of the signals for C₆₀ and C₇₀ were already mentioned in the introduction.^45–49 A comparison between the IR spectrum of the modelled closed cage structures and rugbyballene (the known fullerene, closest by atomic count) is warranted. It is considered that UV radiation in space can break linear and circular carbon chains,⁹¹ hence the form of C₇₂SC without the side chains (C₇₂) is also of interest.

The theoretical IR spectra of C₇₀, C₇₂SC, C₇₀ and C₆₆ are shown in Fig. 8. Table 1 provides frequency comparison for the highest-intensity normal modes of the four structures. Deviations in vibrational frequencies can easily be explained with structural differences: (1) the theoretically assembled nanoparticles contain 4-7- and 8-membered rings, which are not present in rugbyballene, (2) only C₇₀ has order and repeatability in the structural environment of each of its rings and (3) the proposed molecules have “holes” (or simply large rings) ranging from 10 to 13 atoms. For clarification, the smaller, stiffer 4-atomic rings increase the frequency of the normal modes, while the larger, less tense, 7- and 8-membered rings reduce that frequency. Depending on the proportion and arrangement of the size-varying rings in a vibrating molecular fragment, the frequency of the normal mode may shift in either direction.

Fig. 8 Theoretical IR spectrum of the optimized structures: (a) C₇₀, (b) C₇₂SC, (c) C₇₂, and (d) C₆₆. The given frequencies are harmonic.

Table 1 Frequency comparison of high-intensity normal modes of C₇₀, C₇₂SC, C₇₂ and C₆₆

Nanoparticle	Normal mode frequencies [cm^−1]
C70	746	1068	1338	1440	1529	1688
C72SC	684	954	1369	1506	1609	1792
C72	690	1122	1369	1561	1634	1834
C66	703	1000	1412	1514	1560	1694

---

According to general visual comparison, the vibrational spectrum of the theoretically assembled nanoparticles appears similar to that of C₇₀. In terms of frequency shift, C₆₆ comes closest to C₇₀, surprisingly followed by C₇₂SC and finally – C₇₂.

Identification of the highest intensity IR peaks at ∼650–750 cm⁻¹ and ∼1500 cm⁻¹, has been done through the analysis of the computational results. In C₇₀, the discussed normal vibrational modes are located at 746 cm⁻¹ (labelled as V₁) and 1529 cm⁻¹ (labelled as V₂). In the rest of the closed-cage structures, these high intensity vibrations have frequencies close to those of C₇₀. All discussed normal vibrational modes are grouped according to similarity in the atomic displacements, as described in the text below.

We use benzene as an exemplary molecule to study phonons in the vibrations of 6-membered rings. According to displacement vectors, in C₇₀, V₁ has two phonons. One of them is responsible for the out-of-plane component of atomic motion. According to group theory, its symmetry term is B_2g. This energy quantum has a larger contribution to the displacement vectors. The second phonon is of B_2u symmetry and causes the in-plane component of the motion. In each NP, V₁ and the corresponding vibrations affect mostly certain condensed ring fragments, while approximately a half of the structure remains immobile. The symmetry of atomic motion in C₆₆, C₇₂ and C₇₂SC indicates that the corresponding normal modes possess additional phonons.

Atom count parity influences certain characteristics of B_2u ring vibrations. In rings consisting of an even number of atoms, the general direction of the displacement vectors always alternates atom by atom. The difference in the case of an odd number of atoms is that exactly two neighbouring carbons always vibrate in the same general direction. If a bond belongs to rings with different parity, the displacement vectors are hard to predict.

The magnitude of the out-of-plane displacement is similar for V₁ and the corresponding normal modes. The in-plane motion is less pronounced for molecules C₆₆ and C₇₂SC. Visualization of the displacement vectors in the aforementioned vibrations is shown in Fig. 9.

Fig. 9 Displacement vectors for (a) V₁ in C₇₀ and the corresponding harmonic vibrations in (b) C₆₆, (c) C₇₂ and (d) C₇₂SC. The blue atoms have vectors with direction to the inside of the ring and the vectors of the green atoms have the opposite direction.

In C₇₀, V₂ has two phonons. The vibrational symmetry of the first phonon is A_1g, while that of the second is B_2u. Both quanta are responsible for in-plane motion. The second phonon exhibits a larger contribution to the total atomic displacement. In the theoretical closed-cage structures, all normal modes with a frequency of ±300 cm⁻¹ (compared to V₂) contain multiple quanta. As a result, in each harmonic vibration, rarely two rings vibrate in a similar fashion and no rings have pure B_2u or A_1g traits. In all vibrations, most similar to V₂, the dominant phonons are responsible for stretching of certain bonds – not ring vibrations. Atomic displacement is not the only factor used to recognize the normal modes corresponding to V₂. Analysis of the structure of the IR spectrum, including peak position, sequence and intensity, is also employed. Visualization of the displacement vectors in V₂ and its analogical vibrations is shown in Fig. 10.

Fig. 10 Displacement vectors for (a) V₂ in C₇₀ and the corresponding harmonic vibrations in (b) C₆₆, (c) C₇₂ and (d) C₇₂SC. The blue atoms have vectors with direction to the inside of the ring, and the vectors of the green atoms have the opposite direction.

Although the spherical core is not as symmetrical as fullerenes C₆₀ or C₇₀, accounting for the fact that the process initiates with 25 randomly positioned atoms, we consider the final result to be of relevance, in conjunction with the aforementioned reasons.

Considering all factors, even the IR spectra of the modeled, irregular fullerene-like structures do resemble those of C₇₀. Studies involving IR telescopes may detect such improperly shaped molecular systems and misclassify them as established fullerenes.

4. Conclusions

Aggregates of unbound carbon atoms with different initial orientations, densities, and periodicity were used as starting systems to study the formation of carbon clusters and nanophases in space. The stability and evolution of the species were tested using DFTB molecular dynamics. Spheroid carbon nanoparticles resembling imperfect fullerenes were modeled by coupling the resultant clusters through the use of metadynamics.

Isolated carbon aggregations of 20 to 25 atoms form sp²-dominant clusters representing the C-skeleton of PAHs, given a sufficiently high original density. The results are very similar, regardless of initial interatomic orientation. Such clusters exhibit resemblance to fullerene building blocks. If the original density is significantly reduced, similar-sized aggregations result in branched molecules with similar amounts of C atoms in sp and sp² hybridization and no carbons in sp³ state. Such structures have very few rings of low atomic order. If the PBC model has varying cell dimensions, tailored by the cluster size, sp²-dominant multidimensional polymers are formed at higher initial densities. Structural rearrangements are observed during DFTB dynamical simulations. Such simulations tend to show an increase in the quantity of atoms in sp² hybridization and a decrease in the quantity of sp³ hybridized ones. Hence, we also expect such tendency to guide spontaneous transformations in nature. All 0-dimensional forms conform to the experimentally detected types of molecules in space.

The spheroid nanoparticles were assembled with coupling reactions between primordial carbon clusters. One of these clusters was the result of the spontaneous chemistry between carbon atoms randomly arranged in their initial state. Another cluster was formed spontaneously from carbon atoms with a tetrahedral arrangement in the initial system. In each interfragment coupling simulation, most to all C–C binding reactions are barrierless. Structure C₆₆ has the lowest binding energy per atom among the studied molecules. Occasionally, a misfit of atomic positions or a required slight whole-body twist may result in a low barrier for a single coupling reaction. The theoretical IR spectrum of two of the spheroid nanoparticles resembles that of fullerene C₇₀ and it is possible that such imperfect species can be misidentified as known fullerenes in experimental IR studies of carbon signal sources in space.

In this study, we advance a previous hypothesis about the interstellar formation of fullerenes by discovering that graphene is not a required intermediate, and tangled and flat polycyclic structures do not reform only in a direction away from fullerenes. To our knowledge, this is the first research to model a mechanism of formation of fullerene-like structures from primordial, unbound aggregations of carbon atoms.

Data availability

Data for this article are generally included in the article itself as we do not have a separate ESI.

Conflicts of interest

The authors have no conflicts of interest to declare.

Acknowledgements

The authors gratefully acknowledge the financial support from the National Science Fund of Bulgaria under grant KP-06-COST/10 (2023) in the framework of the COST Action CA21126 NanoSpace. This article is based on work from COST Action CA21126 – Carbon molecular nanostructures in space (NanoSpace), supported by COST (European Cooperation in Science and Technology).

References

1. T. J. J. Müller and U. H. F. Bunz, in Functional Organic Materials: Syntheses, Strategies and Applications, ed. T. J. J. Müller and U. H. F. Bunz, Wiley-VCH, 2006.
2. S. Kwok, The synthesis of organic and inorganic compounds in evolved stars, Nature, 2004, 430, 985–991.
3. S. Kwok, Delivery of complex organic compounds from planetary nebulae to the solar system, Int. J. Astrobiol., 2009, 8(3), 161–167.
4. E. F. Van Dishoeck and G. Blake, Chemical evolution of star-forming regions, Annu. Rev. Astron. Astrophys., 1998, 36, 317.
5. S. P. Souza and B. L. Lutz, Detection of C₂ in the interstellar spectrum of Cygnus OB2 Number 12 (IV Cygni Number 12), Astrophys. J. Lett., 1977, 216, L49–L51.
6. P. Ehrenfreund, S. Rasmussen, J. H. Cleaves and L. Chen, Experimentally tracing the key steps in the origin of life: the aromatic world, Astrobiology, 2006, 6(3), 490–520.
7. F. Bertoldi, C. L. Carilli, P. Cox, X. Fan and M. A. Strauss, et al., Dust emission from the most distant quasars, Astron. Astrophys., 2003, 406, L55–L58.
8. J. Qiu, Y. Zhang, J. Zhang and J. Nakashima, Molecules in the Carbon-rich Protoplanetary Nebula CRL 2688, Astrophys. J., Suppl. Ser., 2022, 259, 56.
9. W. Bakke and B. Hrivnak, Searching for Periodicity in 6 Carbon-rich Protoplanetary Nebula, American Astronomical Society Meeting #240, Bull. Am. Astron. Soc., 2022, 54, 6.
10. B. J. Hrivnak and S. Kwok, The Carbon-rich Proto-Planetary Nebula IRAS 22272 + 5435, Astrophys. J., 1991, 371, 631.
11. N. Madhusudhan, M. Agúndez, J. I. Moses and Y. Hu, Exoplanetary Atmospheres – Chemistry, Formation Conditions, and Habitability, Space Sci. Rev., 2016, 205(1–4), 285–348.
12. T. S. Barman, Q. M. Konopacky, B. Macintosh and C. Marois, Simultaneous Detection of Water, Methane, and Carbon Monoxide in the Atmosphere of Exoplanet HR8799b, Astrophys. J., 2015, 804(1), 10.
13. K. B. Stevenson, J. Harrington, S. Nymeyer, N. Madhusudhan, S. Seager and W. C. Bowman, et al., Possible thermochemical disequilibrium in the atmosphere of the exoplanet GJ 436b, Nature, 2010, 464(7292), 1161–1164.
14. M. R. Swain, G. Vasisht and G. Tinetti, The presence of methane in the atmosphere of an extrasolar planet, Nature, 2008, 452(7185), 329–331.
15. B. A. McGuire, R. A. Loomis, A. M. Burkhardt, K. Lee, C. N. Shingledecker and S. B. Charnley, et al., Detection of two interstellar polycyclic aromatic hydrocarbons via spectral matched filtering, Science, 2021, 371(6535), 1265–1269.
16. D. Cruikshank, H. Imanaka, O. Dalle and M. Cristina, Tholins as coloring agents on outer Solar System bodies, Adv. Space Res., 2005, 36, 178–183.
17. F. Raulin, Astrobiology and habitability of Titan, Space Sci. Rev., 2008, 135, 37–48.
18. R. D. Lorenz, K. L. Mitchell, R. L. Kirk, A. G. Hayes, O. Aharonson and H. A. Zebker, et al., Titan's inventory of organic surface materials, Geophys. Res. Lett., 2008, 35, L02206.
19. J. Crovisier, N. Biver, D. Bockelée-Morvan, J. Boissier, P. Colom and C. Dariusz, The chemical diversity of comets: synergies between space exploration and ground-based radio observations, Earth, Moon, Planets, 2009, 105, 267–272.
20. G. Wlodarczak, D. Boucher, R. Bocquet and J. Demaison, The microwave and submillimeter-wave spectrum of methyl chloride, J. Mol. Spectrosc., 1986, 116(1), 251–255.
21. E. C. Fayolle, K. I. Öberg, J. K. Jørgensen, K. Altwegg, H. Calcutt and H. S. P. Müller, et al., Protostellar and cometary detections of organohalogens, Nat. Astron., 2017, 1, 703–708.
22. M. Guélin and J. Cernicharo, Organic Molecules in Interstellar Space: Latest Advances, Front. Astron. Space Sci., 2022,(9), 787567.
23. B. A. McGuire, Census of Interstellar, Circumstellar, Extragalactic, Protoplanetary Disk, and Exoplanetary Molecules, Astrophys. J., Suppl. Ser., 2018, 239(2), 17.
24. P. Swings and L. Rosenfeld, Considerations Regarding Interstellar Molecules, Astrophys. J., 1937, 86, 483–486.
25. A. McKellar, Evidence for the Molecular Origin of Some Hitherto Unidentified Interstellar Lines, Publ. Astron. Soc. Pac., 1940, 52(307), 187.
26. A. A. Penzias, P. M. Solomon, R. W. Wilson and K. B. Jefferts, Interstellar Carbon Monosulfide, Astrophys. J., 1971, 168, L53.
27. S. P. Souza and B. L. Lutz, Detection of C₂ in the interstellar spectrum of Cygnus OB2 Number 12 (IV Cygni Number 12), Astrophys. J. Lett., 1977, 216, L49–L51.
28. H. E. Matthews, W. M. Irvine, P. Friberg, R. D. Brown and P. D. Godfrey, A new interstellar molecule: triearbon monoxide, Nature, 1984, 310(5973), 126.
29. M. Agúndez, J. Cernicharo and M. Guélin, Discovery of interstellar ketenyl (HCCO), a surprisingly abundant radical, Astron. Astrophys., 2015, 577, 6.
30. M. Guelin, S. Green and P. Thaddeus, Detection of the C4H radical toward IRC +10216, Astrophys. J., 1978, 224, L27–L30.
31. P. F. Bernath, K. H. Hinkle and J. J. Keady, Detection of C5 in the circumstellar shell of IRC +10216, Science, 1989, 244, 562–564.
32. J. Cernicharo, C. A. Gottlieb, M. Guelin, T. C. Killian, G. Paubert and P. Thaddeus, et al., Astronomical Detection of H₂CCC, Astrophys. J., Lett., 1991, 368, L39.
33. M. Guelin, N. Neininger and J. Cernicharo, Astronomical detection of the cyanobutadiynyl radical C₅N, Astron. Astrophys., 1998, 335, L1–L4.
34. H. Suzuki, M. Ohishi, N. Kaifu, S. Ishikawa and T. Kasuga, Detection of the interstellar C₆H radical, Publ. Astron. Soc. Jpn., 1986, 38(6), 911–917.
35. L. T. Little, G. H. MacDonald, P. W. Riley and D. N. Matheson, Observations of interstellar HC₅N and HC7N in dark dust clouds, Mon. Not. R. Astron. Soc., 1978, 183, 45–50.
36. J. Cernicharo and M. Guelin, Discovery of the C₈H radical, Astron. Astrophys., 1996, 309, L27–L30.
37. N. W. Broten, T. Oka, L. W. Avery, J. M. MacLeod and H. W. Kroto, The detection of HC₉N in interstellar space, Astrophys. J., Lett., 1978, 223, L105–L107.
38. J. Cernicharo, A. M. Heras, A. G. G. M. Tielens, J. R. Pardo, F. Herpin and M. Guélin, et al., Infrared Space Observatory's Discovery of C₄H₂, C₆H₂, and Benzene in CRL 618, Astrophys. J., 2001, 546(2), L123–L126.
39. B. A. McGuire, A. M. Burkhardt, S. Kalenskii, C. N. Shingledecker, A. J. Remijan and E. Herbst, et al., Detection of the aromatic molecule benzonitrile (c-C₆H₅CN) in the interstellar medium, Science, 2018, 359(6372), 202–205.
40. J. Cernicharo, M. Agúndez, C. Cabezas, B. Tercero, N. Marcelino and J. R. Pardo, et al., Pure hydrocarbon cycles in TMC-1: discovery of ethynyl cyclopropenylidene, cyclopentadiene and indene, Astron. Astrophys., 2021, 649, L15.
41. M. C. McCarthy, K. L. K. Lee, R. A. Loomis, A. M. Burkhardt, C. N. Shingledecker and S. B. Charnley, et al., Interstellar Detection of the Highly Polar Five-Membered Ring Cyanocyclopentadiene, Nat. Astron., 2012, 5, 176–180.
42. J. Cernicharo, M. Agundez, R. I. Kaiser, C. Cabezas, B. Tercero and N. Marcelino, et al., Discovery of two isomers of ethynyl cyclopentadiene in TMC-1: abundances of CCH and CN derivatives of hydrocarbon cycles, Astron. Astrophys., 2021, 655, L1.
43. J. Cernicharo, M. Agundez, R. I. Kaiser, C. Cabezas, B. Tercero and N. Marcelino, et al., Discovery of benzyne, o-C₆H₄, in TMC-1 with the QUIJOTE line survey, Astron. Astrophys., 2021, 652, L9.
44. F. Cataldo, From Elemental Carbon to Complex Macromolecular Networks in Space, in Astrobiology: Future Perspectives, ed. P. Ehrenfreund, L. Becker and J. Blank, Kluwer Academic Publishers, Dordrecht: The Netherlands, 2004. (Reprinted from Astrophysics and Space Science Library, 2004, vol. 305, pp. 97–126).
45. O. Berné, G. Mulas and C. Joblin, Interstellar C₆₀⁺, Astron. Astrophys., 2013, 550, 5.
46. O. Berné, N. L. J. Cox, G. Mulas and C. Joblin, Detection of buckminsterfullerene emission in the diffuse interstellar medium, Astrophys. J., Lett., 2017, 605, L1.
47. J. Cami, J. Bernard-Salas, E. Peeters and S. E. Malek, Detection of C₆₀ and C₇₀ in a Young Planetary Nebula, Science, 2010, 329(5996), 1180.
48. K. Sellgren, M. W. Werner, J. G. Ingalls, J. D. T. Smith and C. Joblin, C₆₀ in Reflection Nebulae, Astrophys. J., Lett., 2010, 722, L54–L75.
49. D. A. García-Hernández, S. Iglesias-Groth, J. A. Acosta-Pulido, A. Manchado, P. García-Lario, L. Stanghellini, E. Villaver, R. A. Shaw and F. Cataldo, The Formation Of Fullerenes: Clues From New C₆₀, C₇₀, and (Possible) Planar C₂₄ Detections In Magellanic Cloud Planetary Nebulae, Astrophys. J., Lett., 2011, 737, L30.
50. Y. Wang, Q. Zeng, X. Du, Y. Gao and B. Yin, The structural, mechanical and electronic properties of novel superhard carbon allotropes: ab initio study, Mater. Today Commun., 2021, 29, 102980.
51. X. Sheng, Q. Yan, F. Ye, Q. Zheng and G. Su, T-Carbon: A Novel Carbon Allotrope, Phys. Rev. Lett., 2011, 106, 155703.
52. D. V. Chachkov and O. V. Mikhailov, Cubic Octa-Carbon: Quantum-Chemical Design of Molecular Structure and Potential Way of Its Synthesis from Cubane, Int. J. Mol. Sci., 2021, 22(21), 12067.
53. T. Gruber and A. Grüneis, Ab initio calculations of carbon and boron nitride allotropes and their structural phase transitions using periodic coupled cluster theory, Phys. Rev. B, 2018, 98, 134108.
54. D. A. Kalchevski, V. Petrov, A. Tadjer and A. Nenov, Impact of hydroxylation on the photophysics of chalcones: insights on the relation between chemical composition and electronic structure, Phys. Chem. Chem. Phys., 2018, 20, 8924–8934.
55. D. Fan, S. Lu, A. A. Golov, A. A. Kabanov and X. Hu, D-carbon: ab initio study of a novel carbon allotrope, J. Chem. Phys., 2018, 149, 114702.
56. D. Li, F. Tian, D. Duan, Z. Zhao, Y. Liu and B. Chu, et al., Modulated T carbon-like carbon allotropes: an ab initio study, RCS Adv., 2014, 4, 17364–17369.
57. S. F. Matar and V. L. Solozhenko, Ultra-hard rhombohedral carbon by crystal chemistry and ab initio investigations, J. Solid State Chem., 2021, 302, 122354.
58. S. Datta, B. Sadhukhan, C. B. Chaudhuri, S. Chakrabarty and A. Mookerjee, Increased metallicity of Carbon nanotubes because of incorporation of extended Stone-Wales’ defects: an ab initio real space approach, Indian J. Phys., 2016, 91, 269–276.
59. G. A. Marchant, M. A. Caro, B. Karasulu and L. B. Pártay, Exploring the configuration space of elemental carbon with empirical and machine learned interatomic potentials, npj Comput. Mater., 2023, 9, 131.
60. G. Pascoli and J. Leclercq, A model for the growth of hydrogenated carbon clusters under interstellar conditions, Astrophys. Space Sci., 1997, 246, 141–152.
61. C. Dubosq, C. Falvo, F. Calvo, M. Rapacioli, P. Parneix and T. Pino, et al., Mapping the structural diversity of C₆₀ carbon clusters and their infrared spectra, Astron. Astrophys., 2019, 625, L11.
62. C. Dubosq, F. Calvo, M. Rapacioli, E. Dartois, T. Pino and C. Falvo, et al., Quantum modeling of the optical spectra of carbon cluster structural families and relation to the interstellar extinction UV bump, Astron. Astrophys., 2020, 634, A62.
63. E. R. Micelotta, A. P. Jones, J. Cami, E. Peeters, J. Bernard-Salas and G. Fanchini, The formation of cosmic fullerenes from arophatic clusters, Astrophys. J., 2012, 761, 35.
64. B. Bhattarai and D. A. Drabold, Amorphous carbon at low densities: an ab initio study, Carbon, 2017, 115, 532–538.
65. B. A. Steele, S. Bastea and I. W. Kuo, Ab initio structural dynamics of pure and nitrogen-containing amorphous carbon, Sci. Rep., 2023, 13, 19657.
66. P. R. Taylor, J. M. L. Martin, J. P. François and R. Gijbels, Ab initio study of the carbon (C₃⁺) cation using multireference methods, J. Phys. Chem., 1991, 95(17), 6530–6534.
67. J. Takahashi, Ab Initio Calculations for Detectable New Isomers of Interstellar Carbon-Chain Radicals C_n(n = 2–8), Publ. Astron. Soc. Jpn., 2000, 52(3), 401–407.
68. T. Chen and A. Li, Synthesizing carbon nanotubes in space, Astron. Astrophys., 2019, 631, A54.
69. A. Chuvilin, U. Kaiser, E. Bichoutskaia, N. A. Besley and A. N. Khlobystov, Direct transformation of graphene to fullerene, Nat. Chem., 2010, 2, 450–453.
70. J. Zhang, F. L. Bowles, D. W. Bearden, W. K. Ray, T. Fuhrer and Y. Ye, et al., A missing link in the transformation from asymmetric to symmetric metallofullerene cages implies a top-down fullerene formation mechanism, Nat. Chem., 2013, 5, 880–885.
71. F. Pietrucci and W. Andreoni, Fate of a Graphene Flake: A New Route toward Fullerenes Disclosed with Ab Initio Simulations, J. Chem. Theory Comput., 2014, 10(3), 913–917.
72. H. W. Kroto and K. McKay, The formation of quasi-icosahedral spiral shell carbon particles, Nature, 1988, 331(6154), 328–331.
73. J. R. Heath, in Fullerenes: Synthesis, properties, and chemistry of large carbon clusters, ed. G. Hammond and V. Kuck, American Chemical Society, 1992.
74. J. M. Hunter, J. L. Fye, E. J. Roskamp and M. F. Jarrold, Annealing Carbon Cluster Ions: A Mechanism for Fullerene Synthesis, J. Phys. Chem., 1994, 98(7), 1810–1818.
75. P. W. Dunk, J. J. Adjizian, N. K. Kaiser, J. P. Quinn, G. T. Blakney and C. P. Ewels, et al., Metallofullerene and fullerene formation from condensing carbon gas under conditions of stellar outflows and implication to stardust, Proc. Natl. Acad. Sci. U. S. A., 2013, 110(45), 18081–18086.
76. C. S. Hansen, E. Peeters, J. Cami and T. W. Schmidt, Open questions on carbon-based molecules in space, Commun. Chem., 2022, 5, A94.
77. J. VandeVondele, M. Krack, F. Mohamed, M. Parrinello, T. Chassaing and J. Hutter, Quickstep: fast and accurate density functional calculations using a mixed Gaussian and plane waves approach, Comput. Phys. Commun., 2005, 167, 103–128.
78. T. D. Kuehne, M. Iannuzzi, M. Del Ben, V. V. Rybkin, P. Seewald and F. Stein, et al., CP2K: an electronic structure and molecular dynamics software package – Quickstep: efficient and accurate electronic structure calculations, J. Chem. Phys., 2020, 152, 194103.
79. M. Elstner, D. Porezag, G. Jungnickel, J. Elsner, M. Haugk and T. Frauenheim, et al., Self-consistent-charge density-functional tight-binding method for simulations of complex materials properties, Phys. Rev. B: Condens. Matter Mater. Phys., 1998, 58(11), 7260–7268.
80. L. Zhechkov, T. Heine, S. Patchkovskii, G. Seifert and H. A. Duarte, An efficient a Posteriori treatment for dispersion interaction in density-functional-based tight binding, J. Chem. Theory Comput., 2005, 1(5), 841–847.
81. J. C. Slater and G. F. Koster, Simplified LCAO method for the Periodic Potential Problem, Phys. Rev., 1954, 94(6), 1498–1524.
82. J. Harris, Simplified method for calculating the energy of weakly interacting fragments, Phys. Rev. B: Condens. Matter Mater. Phys., 1985, 31(4), 1770–1779.
83. M. Elstner, D. Porezag, G. Jungnickel, J. Elsner, M. Haugk and T. Frauenheim, et al., Self-consistent-charge density-functional tight-binding method for simulations of complex materials properties, Phys. Rev. B: Condens. Matter Mater. Phys., 1998, 58, 7260.
84. K. H. Lee, U. Schnupf, B. G. Sumpter and S. Irle, Performance of Density-Functional Tight-Binding in Comparison to Ab Initio and First-Principles Methods for Isomer Geometries and Energies of Glucose Epimers in Vacuo and Solution, ACS Omega, 2018, 3(12), 16899–16915.
85. A. Fihey and D. Jacquemin, Performances of Density Functional Tight-Binding Methods for Describing Ground and Excited State Geometries of Organic Molecules, J. Chem. Theory Comput., 2019, 15(11), 6267–6276.
86. M. Gruden, L. Andjeklović, A. K. Jissy, S. Stepanović, M. Zlatar and Q. Cui, et al., Benchmarking density functional tight binding models for barrier heights and reaction energetics of organic molecules, J. Comput. Chem., 2017, 38(25), 2171–2185.
87. S. Jahangiri, L. Cai and G. H. Peslherbe, Performance of density-functional tight-binding models in describing hydrogen-bonded anionic-water clusters, J. Comput. Chem., 2014, 35(23), 1707–1715.
88. G. Zheng, S. Irle and K. Morokuma, Performance of the DFTB method in comparison to DFT and semiempirical methods for geometries and energies of C₂₀–C₈₆ fullerene isomers, Chem. Phys. Lett., 2005, 412(1–3), 210–216.
89. T. D. Kuhne, M. Krack, F. R. Mohamed and M. Parrinello, Efficient and Accurate Car-Parrinello-like Approach to Born-Oppenheimer Molecular Dynamics, Phys. Rev. Lett., 2007, 98, 066401.
90. A. Laio and M. Parrinello, Escaping free-energy minima, Proc. Natl. Acad. Sci. U. S. A., 2002, 99(20), 12562–12566.
91. A. Laio and F. L. Gervasio, Metadynamics: a method to simulate rare events and reconstruct the free energy in biophysics, chemistry and material science, Rep. Prog. Phys., 2008, 71, 126601.
92. E. Braun, S. M. Moosavi and B. Smit, Anomalous effects of velocity rescaling algorithms: the flying ice cube effect revisited, J. Chem. Theory Comput., 2018, 14, 5262–5272.
93. G. Bussi, D. Donadio and M. Parrinello, Canonical sampling through velocity rescaling, J. Chem. Phys., 2007, 126, 014101.
94. O. Berné and A. G. G. M. Tielens, Formation of buckminsterfullerene (C₆₀) in interstellar space, J. Cryst. Growth, 2011, 109, 401–406.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d4cp02865a

---