Nonlinear optical properties of molecular materials are governed by aggregation and cooperative effects in the condensed phase, requiring first-principles simulations on large molecular assemblies to achieve a realistic description of experiment. While modern electronic-structure theory implementations mitigate key computational challenges through low-memory two-electron integral engines and multifrequency response solvers, practical calculations of third-order nonlinear optical response in extended systems remain prohibitively expensive. The dominant bottleneck arises because conventional approaches require evaluating many independent components of the rank-four hyperpolarizability tensor, each of which is costly, making orientational averaging prohibitively expensive. Here, we introduce an analytic tensor-averaged formulation of cubic response theory within time-dependent Kohn-Sham density functional theory that removes this spatial bottleneck. By performing the orientational average directly at the level of perturbed densities and transformed Fock matrices, the full Cartesian second-order hyperpolarizability tensor is never constructed. Instead, the isotropic second-order hyperpolarizability is obtained directly. Exploiting the linearity of the Fock-matrix construction, this approach eliminates redundant spatial components and reduces the number of required exchange-correlation kernel integrations in the electronic quartic- and cubic-Hessian contractions by 72 and 90%, respectively. The resulting methodology enables cubic-response simulations of third-harmonic generation in large molecular aggregates with exchange-correlation functionals spanning multiple rungs of Jacob's ladder. Applications to oligo(thiophene-benzothiadiazole) (OTBP) clusters reveal pronounced aggregation-induced dampening of the third-harmonic response and demonstrate scalability to systems with nearly 5800 contracted basis functions. The present framework provides a quantitative tool for disentangling intrinsic molecular effects from supramolecular and morphological contributions, and for assessing the potential for further performance gains through crystal-structure engineering. In this way, theory can directly guide experimental efforts by identifying packing motifs that maximize cooperative enhancement and those that are detrimental. By resolving the missing spatial-domain optimization in the cubic response theory, the present tensor-averaged formulation shifts the dominant computational cost back to the frequency-dependent response solver, enabling routine material-level modeling of cooperative third-order nonlinear optical phenomena.

A software implementation of analytic geometric derivatives of electron-repulsion integrals up to second order is presented for the modeling of vibrational spectroscopies at the level of first-principles Kohn–Sham density functional theory (DFT). In line with the general goals of the VeloxChem program, it targets efficient execution in high-performance computing environments with a hybrid MPI/OpenMP parallelization model and is based on the technique of automatic C++ code generation for high versatility. Gradient calculations scale identically with conventional Fock matrix constructions, and also with the prefactor taken into account, the computational cost of the gradient is significantly lower than that of the self-consistent field (SCF) optimization of the reference state. The Hessian calculation shows a scaling of N3.5 with N being the number of contracted Gaussian basis functions. The computational bottleneck in the Hessian calculation is the solving of the coupled-perturbed Kohn–Sham equations that with VeloxChem can be offloaded to GPU-accelerated nodes. The large-scale virtues of the presented software module are demonstrated by the DFT/B3LYP calculation of the IR spectrum of the entire ubiquitin protein with 1,152 atoms in the quantum mechanical (QM) region and TIP3P water in the molecular mechanics (MM) region. The simulated amide I band shows to be in excellent agreement with experiment.

The calculations of electronic g-tensors, one of the most important parameters in electron paramagnetic resonance spectroscopy, are carried out for dangling bonds ( DBs ) introduced into hydroxylated and aminated diamond nanoparticles, or nanodiamonds ( NDs ) , of different shapes and sizes. Regarding the shapes of NDs, octahedral, cubic, and tetrahedral model systems are used, while the impact of the change in size is inspected by increasing octahedral ND from C 35 to C84. The results for single DBs reveal that tetrahedral NDs exhibit the widest variation range of the isotropic g-shift values for both surface functionalization schemes, whereas the isotropic g-shifts of octahedral and cubic NDs tend to strongly overlap. On the other hand, if one treats NDs as an ensemble of nanoparticles constituting a sample, the isotropic g-shifts arithmetically averaged over all available DBs show that tetrahedral ND with hydroxylated surface possesses a significantly higher value than the rest of the considered systems. However, applying the Boltzmann distribution results in a substantially lower value for cubic ND. In contrast, aminated NDs do not demonstrate average values that stand out from the others, irrespective of the analysis method employed. Overall, in addition to the comprehensive magnetic properties, the obtained data also provide interesting details on the formation of DBs in hydroxylated and aminated NDs.

Computer simulations play an essential role in the interpretation of experimental multiphoton absorption spectra. In addition, models derived from theory allow for the establishment of "structure-property" relationships. This work contributes to these efforts and presents the results of an analysis of two- and three-photon absorptions for a set comprising 450 conjugated molecules performed at the CAM-B3LYP/aug-cc-pVDZ level. The molecular set is composed of organoboron dyes presenting various core topologies combined with a palette of conjugated linkers giving donor-acceptor architectures. The charge-transfer character of the investigated structures is manifested by the presence of the low-lying electronic excited state. The multiphoton excitation to the state in question is intense and significant from an application point of view. The analysis performed in this work clearly demonstrates that there is a strong correlation between the intensities of the two- and three-photon transitions to the lowest intramolecular charge-transfer state, hinting that developed design rules aiming at maximizing two-photon absorption efficiency will also be useful in designing three-photon absorbers. As part of this study, we also performed two-photon absorption calculations using the coupled-cluster RI-CC2 model with the aug-cc-pVDZ basis set for 450 molecules to guide the selection of the density functional approximation.

In this Letter, we present a pioneering analysis of the density functional approximations (DFAs) beyond the generalized gradient approximation (GGA) for predicting two-photon absorption (2PA) strengths of a set of push-pull π-conjugated molecules. In more detail, we have employed a variety of meta-generalized gradient approximation (meta-GGA) functionals, including SCAN, MN15, and M06-2X, to assess their accuracy in describing the 2PA properties of a chosen set of 48 organic molecules. Analytic quadratic response theory is employed for these functionals, and their performance is compared against the previously studied DFAs and reference data obtained at the coupled-cluster CC2 level combined with the resolution-of-identity approximation (RI-CC2). A detailed analysis of the meta-GGA functional performance is provided, demonstrating that they improve upon their predecessors in capturing the key electronic features of the π-conjugated two-photon absorbers. In particular, the Minnesota functional MN15 shows very promising results as it delivers pleasingly accurate chemical rankings for two-photon transition strengths and excited-state dipole moments.

The eChem project features an e-book published as a web page (10.30746/978-91-988114-0-7), collecting a repository of Jupyter notebooks developed for the dual purpose of explaining and exploring the theory underlying computational chemistry in a highly interactive manner as well as providing a tutorial-based presentation of the complex workflows needed to simulate embedded molecular systems of real biochemical and/or technical interest. For students ranging from beginners to advanced users, the eChem book is well suited for self-directed learning, but workshops led by experienced instructors and targeting student bodies with specific needs and interests can readily be formed from its components. This has been done by using eChem as the base for a workshop directed toward graduate students learning the theory and practices of quantum chemistry, resulting in very positive assessment of the interactive nature of this framework. The members of the eChem team are engaged in both education and research, and as a mirroring activity, we develop the open-source software upon which this e-book is predominantly based. The overarching vision and goal of our work is to provide a science- and education-enabling software platform for quantum molecular modeling on contemporary and future high-performance computing systems, and to document the resulting development and workflows in the eChem book.

The geometry optimization of 30 paramagnetic defects, including biomedically attractive nitrogen-, silicon-, germanium-, and nickel-related color centers, is performed after their incorporation into hydrogenated nano -diamond (ND) of C84H64 size. The main aim is to examine the effectiveness of the low-cost methods, namely, PBEh-3c, r2SCAN-3c, B97-3c, HF-3c, and GFN2-xTB, in reproducing the geometries of these defects basing on the similarity between the results of the subsequent electronic g-tensor calculations. It is revealed that the overall performance of PBEh-3c, r2SCAN-3c, and B97-3c is very alike and can be considered as good, however, none of these "3c" approaches is able to cope with all tested geometries. The results of HF-3c, on the other hand, are disappointing, as this method is outperformed by computationally much more lighter GFN2-xTB. Additional calculations carried out for dangling bonds introduced into hydroxylated and aminated NDs show that all low-cost methods perform reasonably well for this type of defect but the largest quantitative discrepancies once again are demonstrated by HF-3c. The obtained findings lay the foundations for the future studies of larger NDs with the purpose to figure out the magnetic properties dependence on the size of NDs or defect positions within NDs.

This study employed a coarse-grained Monte Carlo (MC) simulation to investigate the radiation-induced polymerisation of methacrylic acid (MAA) in an aqueous solution. This method provides an alternative to traditional kinetic models, enabling a detailed examination of the micro-structure and growth patterns of MAA polymers, which are often not captured in other approaches. In this work, we generated multiple clones of a simulation box, each containing a specific chemical composition. In these simulations, every coarse-grained (CG) bead represents an entire monomer. The growth function, defined by the chemical behaviour of interacting substances, was determined through repeated random sampling. This approach allowed us to simulate the complex process of radiation-induced polymerisation, enhancing our understanding of the formation of poly(methacrylic acid) hydrogels at a microscopic level; while Monte Carlo simulations have been applied in various contexts of polymerisation, this study’s specific approach to modelling the radiation-induced polymerisation of MAA in an aqueous environment, utilising the data obtained by quantum chemistry modelling, with an emphasis on micro-structural growth, has not been extensively explored in existing studies. This understanding is important for advancing the synthesis of these hydrogels, which have potential applications in diverse fields such as materials science and medicine.

We present a selected set of exemplifying applications of the novel polarizable coarse-grained model [see the first part] to various outstanding problems in the physics and chemistry of nanoparticles: electrostatic potential around silver and gold nanoparticles; spontaneous and external electric field-driven self-organization of gold and silver nanoparticle systems; and physisorption of carbon dioxide on titanium dioxide nanoparticles decorated with a gold catalyst. In the first application, the developed model has shown capabilities of predicting long-range potential with accuracy comparable to the tight-binding density functional theory methods. Furthermore, the electrostatic potential analysis in hot spot regions allowed us to identify twin defect lines in a silver nanostar as a promising candidate for an enhancer in surface-enhanced Raman spectroscopy. In the second application, the developed model has facilitated the elucidation of the microscopic mechanisms responsible for the self-organization of gold and silver nanoparticles. Analysis of Monte-Carlo simulations established that the self-organization process is driven by van der Waals interactions in the absence of an external electric field, and that it becomes gradually driven by electrostatic interactions in the presence of an external electric field with increasing strength of the external electric field. In the third application, the developed model combined with Monte-Carlo simulations has identified the dominant mechanism responsible for carbon dioxide transfer to the catalytic sites. Analysis of the obtained results indicates that surface diffusion is the dominant mechanism for the transport of carbon dioxide to the catalytic sites, and only in exceptional situations, direct physisorption becomes a competitive mechanism with the surface diffusion mechanism. These successful applications of the developed model indicate its wide range of applicability to various problems in the chemistry and physics of nanoparticles.

We present a polarizable coarse-grained model for metal, metal oxide, and composite metal/metal oxide nanoparticles with well-defined crystalline surfaces. The developed model uses a low-resolution polarizable “surface beads” representation of the nanoparticle's geometry and pairwise cross nanoparticle potential consisting of van der Waals and electrostatic interaction terms. The electrostatic interaction term of the cross nanoparticle potential incorporates a crucial physical aspect of electrostatic interaction into the metal and metal oxide systems, such as induced surface charges, making it possible to explore the nanoparticles’ behavior in complex environments as well as investigate the interplay between electrostatic and van der Waals interactions in nanoparticle systems. The iterative stability, computational scaling, and performance of the presented model was tested on selected systems of gold, titanium dioxide, and composite gold/titanium dioxide nanoparticle systems. The model exhibits robust iterative stability and is able to converge the charge equilibration equation for fluctuating induced charges and dipoles within 10-60 “tug-tow” iterations in challenging situations, like crowded nanoparticle systems or nanoparticle systems in extreme external electric fields. The computation scaling of the presented model is semi-linear with respect to the number of nanoparticles in the system. It slightly varies depending on the size distribution of nanoparticles in a specific nanoparticle system. The computation cost of the model is significantly lower than that of conventional atomistic polarizable force field models and enables the treatment of large nanoparticle systems that are beyond the reach of currently existing atomistic force field models.