Dalton Project: A Python platform for molecular- and electronic-structure simulations of complex systems

The Dalton Project provides a uniform platform access to the underlying full-ﬂedged quantum chemistry codes Dalton and LSDalton as well as the PyFraME package for automatized fragmentation and parameterization of complex molecular environments. The platform is written in Python and deﬁnes a means for library communication and interaction. Intermediate data such as integrals are exposed to the platform and made accessible to the user in the form of NumPy arrays and result data are extracted, analyzed, and visualized. Complex computational protocols that may for instance arise due a need for environment fragmentation and conﬁguration-space sampling of biochemical systems are readily assisted by the platform. The platform is designed to host additional software libraries and will serve as a hub for future modular software development eﬀorts in the distributed Dalton community.


Denmark
(Dated: 26 March 2020) The Dalton Project provides a uniform platform access to the underlying full-fledged quantum chemistry codes Dalton and LSDalton as well as the PyFraME package for automatized fragmentation and parameterization of complex molecular environments.
The platform is written in Python and defines a means for library communication and interaction. Intermediate data such as integrals are exposed to the platform and made accessible to the user in the form of NumPy arrays and result data are extracted, analyzed, and visualized. Complex computational protocols that may for instance arise due a need for environment fragmentation and configuration-space sampling of biochemical systems are readily assisted by the platform. The platform is designed to host additional software libraries and will serve as a hub for future modular software development efforts in the distributed Dalton community.

I. INTRODUCTION
More than 20 years have passed since the first version of the Dalton program 1 was released as the result of merging the separate HERMIT, SIRIUS, ABACUS, and RESPONS codes that implemented one-and two-electron integrals, wave functions, energy derivatives, and response theory, respectively. Later, the adopted monolithic code development structure turned out prohibitively difficult to sustain and it was interrupted with the release of the atomic-orbital (AO) based linear-scaling initiative as a separate executable named License version 2.1 (LGPLv2.1). In the present work, we will briefly recapitulate the functionalities of the codes and detail some of the developments provided in Dalton suite releases from 2015 until today-including the Dalton2020 release. With inspiration from the Molecular Sciences Software Institute (MolSSI) project, 4,5 we also take the opportunity to initiate a transition in the Dalton software engineering practices and we signal this paradigm change by referring to the Dalton community effort as the Dalton Project (DP) initiative. 6 From the developer's perspective, we are taking steps to make it easier to develop, sustain, and maintain a large general-purpose software ecosystem for first-principles quantum molecular modeling of complex systems and, from the user's perspective, we are modifying the design of the user interface to enable new access and interaction patterns.
The general design strategy for the DP platform is that of software modularity 7-9 and based on a hybrid programming language approach as illustrated in Fig. 1  Python bindings, e.g., through CFFI (C Foreign Function Interface) 13 or pybind11, 14 or direct Python module import. In this scheme, we view the Dalton and LSDalton executables as libraries serving the DP platform and, although further modular library decomposition would be desirable, it is hampered by code legacy and entanglement. More important than offering this new perspective, however, the DP platform encourages future development to be made in the form of modules with clear and specific tasks (or subtasks) and which Within the field of quantum chemistry, the adoption of more modern software engineering strategies with APIs written in Python is in vogue at the moment, and we have been strongly influenced by (i) the Psi4NumPy project that exposes efficient computational kernels from the Psi4 program 15 to enable quick NumPy prototyping of novel science 16 and (ii) the PySCF program that, primarily in Python, implements self-consistent field (SCF) and post-Hartree-Fock (post-HF) electronic structure theory for finite and periodic systems. 17 Moreover, a source of inspiration as well as practical experience for the present work is provided by the VeloxChem project (and program) 18 that, with a hybrid Python/C++ programming model, implements real and complex response theory 19 at the SCF level of theory for execution in high-performance computing (HPC) cluster environments. In VeloxChem, Python is used for a split message passing interface (MPI) communicator management of large-scale distributed hardware resources with an anticipation of heterogeneous cluster nodes to become a future reality. Without noticeable sacrifice in computational efficiency or program execution stability, the higher level quantum chemical methods and iterative linear response equation solvers are implemented in Python with use of NumPy and underlying threaded math kernel libraries. With this as background, we have gained sufficient confidence to steer our project into a new direction as far as software engineering practices are concerned.
Our presentation is organized as follows: In Section II, we briefly mention some of the key features in Dalton and LSDalton that have already been presented 3 and provide a more detailed description of novel functionalities that have been added thereafter. Moreover, we also present the features of the PyFraME package 20 for the handling of complex chromophore environments. In Section III, we give a more comprehensive view of the design of the DP platform as well as provide a concrete illustration of new library-access patterns and program execution practices for the user. In Section IV, we present six concrete examples of DP platform runs before closing with an outlook into the future for the Dalton Project.

II. CAPABILITIES OF DP PLATFORM LIBRARIES
A. Dalton and LSDalton up until 2014 In 2014, a presentation of the Dalton program system, including the Dalton and LSDalton codes, was published 3 and functionalities listed in this presentation are of course still 5 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.
available and therefore only briefly mentioned here. The two codes are primarily written in Fortran but parts involving density functional theory (DFT) kernels are mostly written in C.
In Dalton, the routines for correlated wave-function calculations are implemented only for serial execution-but can be linked to standard threaded linear algebra libraries-whereas the self-consistent field (SCF), i.e., HF and DFT, routines are implemented for parallel execution using MPI. LSDalton, on the other hand, comes with a native hybrid OpenMP/MPI parallelization scheme that enables shared memory data handling on central processing unit (CPU) sockets and/or compute nodes. None of the two codes come with support for hardware acceleration, such as general-purpose graphical processing units (GPUs).
The common foundation for the Dalton and LSDalton quantum chemistry programs is that of a nonrelativistic Hamiltonian, basis sets expanded in localized Gaussian AOs, and multi-electron reference states expanded in spin-restricted determinants or configurationstate functions. Relativistic corrections to the zeroth-order one-electron Hamiltonian are available in Dalton in terms of the spin-free second-order Douglas-Kroll-Hess (DKH2) Hamiltonian and effective-core potentials (ECPs). As a perturbative correction to the Hamiltonian, Dalton also offers an implementation of the full Breit-Pauli spin-orbit operator.
LSDalton provides efficient acceleration techniques for SCF-based property calculations and an implementation of the linear-scaling divide-expand-consolidate (DEC) scheme for second-order Møller-Plesset (MP2) and coupled cluster (CC) energy calculations. The code was initially developed to alleviate the restrictions of the Dalton code for calculations on large systems by introducing linear-scaling AO-based SCF and response capabilities based on an exponential ansatz of the AO density matrix.
Dalton provides implementations of most of the standard electronic-structure methods, including SCF, MP2, a hierarchy of CC methods (CC2, CCSD, CCSDR(3), CC3, and CCSD(T)), configuration interaction (CI), and multi-configurational SCF (MCSCF) based on the generalized active space (GAS) concept. MCSCF wave functions are optimized with a robust trust-region-based second-order approach.
Molecular gradients and Hessians are determined analytically for SCF and MCSCF reference states and analytic gradients are also available at the levels of MP2, CC2, CCSD, and CCSD(T). In the absence of analytic gradients and Hessians, Dalton can determine these quantities by numerical differentiation and thereby offers an extensive functionality for exploring potential energy surfaces. The combination of geometric and electric-field per-6 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5144298
turbations allows for calculations of infrared (IR) and Raman intensities. Analytic linear and nonlinear response functions describing the interactions with external and internal (in general time-dependent) electromagnetic fields are implemented for the entire selection of electronic-structure methods and enable simulations of a plethora of spectroscopies, too rich to be listed here. At this time, Dalton also included the means for structure-less and atomistic descriptions of chromophore environments through the polarizable continuum model (PCM) and polarizable embedding (PE) model, respectively.

B. Added features in Dalton
The functionalities of Dalton have been expanded in several directions. Here we provide a summary of selected new features. To bring some structure and order into these developments, we have chosen to divide them into the three categories: (i) electronic-structure theory, (ii) spectroscopy simulations, and (iii) environment modeling. In the first one, we list general quantum-chemical method developments providing new means to describe the electronic structure of ground and excited states. In the second, we describe developments that are more specifically targeting and enabling simulations of certain spectroscopies. Such simulations are connected with certain electronic-structure theory methods and typically also environment models, but the primary objective of the development has been the spectroscopy at hand. In the third, we present approaches aimed at improving the effective description of the chromophore environment. These developments are of course made in combination with specific electronic-structure methods but the environment is at focus.

Electronic-structure theory
Based on a range-separated Hamiltonian as proposed by Savin, 21,22 a rigorous combination can be made of wave-function and density-functional theories for the treatment of the wave function at the same time as dynamic electron correlation can be effectively accounted for by means of DFT with its low computational cost. The resulting method is referred to as multi-configurational short-range DFT (MC-srDFT) and it is available for closed-shell and open-shell systems. 28 Apart from calculations of energies, linear-response properties are available for both singlet and triplet perturbations. [30][31][32] More details are provided in Section IV B where an example is provided in terms of the calculation of the ultraviolet-visible (UV/Vis) absorption spectrum of a retinylidene Schiff base chromophore.
Using Löwdin's inner projection in conjunction with a one-and two-electron excitation operator manifold and an MP2 reference state, the second-order polarization propagator approximation (SOPPA) arises as a means to address the electronic structure of excited states. Going instead toward approximations of the SOPPA model, the random phase approximation with second-order non-iterative doubles corrections model (RPA(D)) 37 has been extended to triplet excitations 38 and a similarly derived higher RPA with non-iterative doubles corrections model (HRPA(D)) has been implemented. 38 The RPA(D) and HRPA(D) models are enabled for calculations of not only transition properties but also linear response functions. 39 Furthermore, in regard with CC approaches, Dalton also offers a new and more efficient implementation of the CC3 model for ground-and excited-state energies 40 -it does not support full Abelian point-group symmetry but reduces the computational cost.

Spectroscopy simulations
Applying the Liouville equation to pure states in the density-matrix formalism of quantum mechanics has been shown equivalent to applying the Ehrenfest theorem to state-transfer operators in Hilbert space and thereby leading to a means to phenomenologically introduce relaxation mechanisms into wave-function theories. 19,41,42 The resulting complex polarization propagator (CPP) theory defines frequency-dependent response functions for exact and approximate states that are physically sound in all regions of the spectrum, resonant as well as conventional nonresonant, and also X-ray as well as conventional UV/Vis. These re-8 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. Core excitation processes are associated with large valence-electron relaxation and polarization effects that in a polarization propagator or response theory approach require multi-electron excited configurations to be properly accounted for. 53 Along this line, Dalton provides a hierarchy of CC methods to model a variety of X-ray spectroscopies including near-edge X-ray absorption fine structure (NEXAFS), 54-57 photo-electron spectroscopy (PES), 56-59 transient X-ray absorption spectroscopy (TRXAS), 60,61 and resonant inelastic Xray scattering (RIXS). 62 The referred-to hierarchy of CC methods includes the CCS, CC2, and CCSD levels of theory, but core-excitation and core-ionization energies are also available for the CCSDR(3) and CC3 approximations. Both singlet and triplet excited states are encompassed and the latter are obtained in a spin-adapted formalism. 61 The core-valence separation (CVS) approximation has been made available to decouple core and valence excited states. It can be applied either at the excited-state level only 56,57 or both during the determination of the ground state and excited states in a frozen-core variant (fc-CVS). 63 An example illustration of a DP platform XAS calculation using the CVS approximation is provided in Section IV F.

Modeling of chromophore environments
The capabilities in Dalton for including effects from a molecular environment have been extended in several directions. The PCM for efficient modeling of bulk solvent effects can now also be performed at the SOPPA level. 64  This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.
with the MC-srDFT method.
The PE model 65,66 is a fragment-based quantum-classical scheme designed for efficient and accurate inclusion of environment effects in calculations of spectroscopic properties of large and complex molecular systems. The environment is included effectively through an embedding potential whose parameters consist of distributed multipoles and polarizabilities, both of which are derived from quantum-mechanical calculations on the individual fragments that make up the environment. The PyFraME package, which is made available on the DP platform and is described in Section II D, can be used to automatize the generation of the embedding-potential parameters. The PE model is implemented in the Polarizable Embedding library (PElib) 67 based on an AO density-matrix-driven formulation, which facilitates a loose-coupling modular implementation in host programs. The PElib was included in the Dalton2013 release but at that time it was limited to PE-HF and PE-DFT. 65,66 Since then the implementation has been extended to PE-CC (specifically PE-CC2, PE-CCSD, and PE-CCSDR(3)), 68 PE-MCSCF, 69 PE-MC-srDFT, 70 and PE-SOPPA. 71 The Dalton2020 release supports linear-, quadratic-, and cubic-response properties for PE-HF/DFT, 66 while PE-CC is limited to linear-and quadratic-response properties, and only linear-response properties are available for PE-MCSCF and PE-MC-srDFT. For PE-HF/DFT it is also possible to compute properties based on resonant-convergent response. 72 London AOs (LAOs) are supported for magnetic linear-response properties that involve a single derivative with respect to a magnetic field. 73 The capabilities have also been extended to enable analytic quantummechanical molecular gradients at the PE-HF/DFT level, thus enabling geometry optimization of the core quantum region embedded in a fixed polarizable environment. 74 Local-field effects may also be included in PE-HF/DFT calculations where they are termed effective external field (EEF) effects. 75 The PE model, as with other quantum-classical embedding models, does not include Pauli repulsion between the chromophore and its environment. Such models can therefore suffer from so-called electron spill-out, where the electron density of the chromophore leaks This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5144298
out into the environment, thus causing an over-stabilization of the ground state and, in particular, the excited states of the embedded chromophore. Negatively charged chromophores or excited states of even partial Rydberg-like character are especially susceptible. 83,84 The polarizable density embedding (PDE) model has been formulated to improve the electrostatic interactions between the chromophore and its environment, and to address the electron spill-out issue. 85,86 In this model, the permanent charge distribution of the fragments in the environment are described by their full electronic densities, thus avoiding divergences of the multipole expansions, while still keeping the distributed polarizabilities to efficiently account for polarization effects. In addition, the PDE model contains a Huzinaga-Cantu-like projection operator 87 that models Pauli repulsion effects and thereby effectively prevents electron spill-out. The PDE model has been implemented in PElib using the same AO density-matrix-driven formulation as for the PE model. The PDE model can therefore straightforwardly be combined with the same DFT and wave-function methods as the PE model both in terms of ground-state and response calculations, with the exception of LAOs and analytic molecular gradients.
The latest Dalton release has also received basic frozen density embedding (FDE) 88,89 capability. The FDE implementation enables import of a static embedding potential that has been pre-calculated on a numerical integration grid by another code that implements FDE. 90,91 A matrix representation of the embedding potential is constructed based on the grid and added to the one-electron Fock matrix. The implementation can thus be used with all available DFT and wave function methods in Dalton. 92,93 C. Added features in LSDalton

Integral evaluation
Integrals sit at the heart of any quantum chemistry program, both when it comes to computational performance and available methods and properties. The development of an efficient and flexible integral-evaluation code has therefore been essential to the development of LSDalton. Since 2014, 3 four main integral developments have been added: high-order derivative integrals (HODI), integrals and differentiated integrals for embedding techniques involving interaction with point charges and higher-order multipoles, acceleration of the exchange 11 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.
contribution through developments of the auxiliary-density-matrix method (ADMM), 94 and interface with the XCFun library of DFT exchange-correlation (XC) functionals. 95,96 The XCFun library is based on forward-mode automatic differentiation 97 and can therefore generate arbitrary-order derivatives of these functionals. for the projection to the smaller basis and GGA correction functional options, 99 and with tailored auxiliary basis sets (admm-n) for the pcseg-n and aug-pcseg-n basis sets. 100

Exploiting the locality of electron correlation
The Divide-Expand-Consolidate (DEC) 101-103 strategy employs highly local orbitals 104,105 to recover the inherent locality of dynamical correlation for large molecules in a linearscaling fashion. Over the last years, the DEC framework has been extensively developed, and now includes resolution-of-the-identity (RI) accelerated MP2 (DEC-RI- This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5144298
through the local framework for calculating excitation energies (LoFEx) [112][113][114] and the correlated natural transition orbital framework for low-scaling excitation energy (CorNFLEx) approach. 115 Due to the embarrassingly parallel nature of the DEC scheme, highly parallel scaling is possible. As an example, a cc-pVDZ DEC-RI-MP2 gradient calculation of the insulin molecule (787 atoms and 7,604 basis functions) finished in less than 10 hours using 32,000 cores (2,000 nodes, each with 16 cores, of the Titan supercomputer). 107

Molecular properties
Several property developments have been undertaken since 2014, including quasi-Newton transition-state optimization, the high-order path-expansion (HOPE) 116 method for improved geometry-optimization steps, automated counterpoise correction, and the samenumber-of-optimized-parameters (SNOOP) 117,118 scheme as an improved alternative to the counterpoise correction, and nuclear-selected NMR shielding, 119 to mention a few.
On a longer-term development line, LSDalton has been interfaced with OpenRSP, 120 to allow, in principle, arbitrary-order molecular properties.
OpenRSP is an open-ended response-theory library that that manages the generation and solution of the response equations needed for the evaluation of arbitrary-order response properties. Currently, the interface to LSDalton enables calculation of (mixed) general-order electrical and geometrical properties for HF and up to fourth order for DFT through an interface to the XCFun and XCint 121 libraries, including, for example, IR, Raman, and hyper-Raman spectroscopic properties, and molecular gradients, Hessians, cubic, and quartic forces. The capabilities of the LSDalton/OpenRSP/XCint/XCFun combination is illustrated for the calculation of IR and Raman spectra of benzene through the DP platform in Section IV D.
For modeling of solvent effects, the PCMSolver 122 library for continuum electrostatic solvation has been interfaced to LSDalton. This implementation is available for SCF and electric-dipole response properties up to fourth order.

D. PyFraME
PyFraME 20 is a Python package providing a framework for managing fragment-based multiscale embedding calculations. The basic principle of embedding models in quantum 13 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. chemistry is the division of a molecular system into two domains: a central core region and its environment. The core region is described at the highest level of theory using DFT or a wave function method while the effects from the environment are included effectively through an embedding potential. To manually set up embedding calculations of large and complex molecular systems can be highly complicated, tedious, and error-prone. This is especially true when considering that configuration-space sampling, e.g., through molecular dynamics (MD) simulations, is usually required and which in turn means that the procedure has to be repeated many times.
The highly flexible PyFraME package automatizes workflows, starting from initial molecular structure to embedding potential, enabling the user to easily setup a multilayered description of the environment. Each layer can be described either by a standard embedding potential, i.e., using a predefined set of parameters, or by deriving the embedding-potential parameters based on first-principles calculations. For the latter, a fragmentation method is used to subdivide large molecular structures into smaller computationally manageable fragments. The number of layers, as well as the composition and level of theory used for each layer, can be fully customized.
The basic workflow consists of three main steps. First, a molecular structure is given as input. Currently, PyFraME supports input files in the PDB format. The input file reader extracts information about the structure and composition of the system and it also defines the basic units of the system, i.e., fragments. Small molecules would typically constitute a fragment on their own but larger molecules are usually broken down into small computationally manageable fragments. For example, for proteins a fragment would usually consist of an amino-acid residue while for nucleic acids it could be a nucleotide. The molecular system to be used for the embedding calculation is then built by extracting subsets from the full list of fragments according to specified criteria, such as name, chain ID, distance, or a combination thereof, and placed into separate regions. As mentioned above, any number of regions may be added and each can be fully customized. Once the system has been built, the final step is the derivation of the embedding potential. Depending on the specifics, it may involve a large number of separate calculations on the individual fragments in order to compute the embedding-potential parameters. For large molecules, where the parameters can not be computed directly, PyFraME uses a fragmentation method based on the molecular fractionation with conjugate caps (MFCC) approach 123 to derive the parameters. The

14
This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

III. DP PLATFORM DESIGN AND FEATURES
The ultimate goal of the DP platform is to establish a flexible, robust, and uniformly accessible environment that can be used for both large-scale applications and to facilitate the development of novel methodology. The challenges for quantum chemistry today are far more complex than earlier, both concerning the complexity of the chemical systems, adaptation to HPC facilities, and the number of tools and approaches needed for applications. As a result, it becomes essential to be able to easily combine the tools and approaches in meaningful ways. The main motivation of the Dalton Project is to provide a platform that can be used to combine the functionality of the tools and methods that are developed by the individual research groups in the Dalton community.
For a long time, we have relied on a monolithic codebase (first Dalton and later also LSDalton) for the development of new computational methodology. These programs have served us well in the past and we expect this to continue into the foreseeable future. It is clear, however, that the codebase has accumulated substantial technical debt. The tight coupling between the software modules is particularly problematic because it complicates optimization and modernization of even small pieces of code. Moreover, implementation of new methodology often requires unnecessarily high efforts and easily leads to additional technical debt. The risks of relying on a monolithic codebase are especially high when the codebase is maintained by a scientific community such as ours whose primary goal is to perform research. In recognition of this, and the fact that individual groups have different research aims and preferences in terms of software development, we have in later years moved towards a more distributed codebase. This has resulted in the development of a series of software libraries, such as Gen1Int, 126,127 OpenRSP, 120 PCMSolver, 122 PElib, 67 QcMatrix, 128 XCFun, 95,96 and XCint. 121 This has, to some degree, alleviated the problem of the monolithic codebase for some developments but the main issues remained.
We have now taken the next step and moved completely to a distributed codebase with the Dalton Project, whose main task is to integrate and provide interoperability between the individual software libraries that are developed and maintained by the different nodes in our community. At the same time, however, we acknowledge that there is vast functionality developed in our community during the last decades that we do not wish to abandon, and which is primarily implemented in Dalton and LSDalton. The DP platform thus has to accommodate a wide variety of software, from large monolithic programs with a wide range of features to small libraries that provide very specific functionality. The design of the platform has to take this into account in a sustainable manner, so that it can act as a platform for present-day use cases and, importantly, for future developments based on modern software engineering practices. Moreover, the DP platform must be able to exploit current HPC facilities and be prepared for the upcoming exascale supercomputers.
To meet our goals and requirements, we devised a platform structure that is illustrated in Fig. 1. At the top, we have the DP platform itself, written in pure Python (3.6+), that interfaces to external libraries through different communication mechanisms. Python was chosen as the platform language because of its extensibility, emphasis on code readability, and comprehensive standard library, as well as the large number of specialized libraries. We note here that the term library is used liberally to signify any type of software that can be LSDalton but we will gradually move towards using loosely coupled modules written in pure Python or hybrid Python and Fortran, C, or C++, with the hybrid approach being used for the more compute-intensive numerical tasks. Initially, the DP platform will interface Dalton and LSDalton as well as PyFraME, all of which have been described in the previous section. In the immediate future, we expect that many of the aforementioned libraries that are currently interfaced to Dalton and LSDalton will be interfaced directly to the DP platform.
The DP platform is a Python package with an API that consists of a set of classes and functions used to setup molecular systems, perform numerical calculations, and process data. Usage of the platform would typically consist of three stages used in succession, namely 16 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. setup, compute, and data processing. By separating the compute and data processing stage, the DP platform may be easily employed in large-scale application workflows, in which a number of calculations are typically run first, by submitting them to a queuing system on a supercomputer, and subsequently data manipulation, analysis, and visualization is performed.

PLEASE CITE THIS ARTICLE AS
The setup stage consists of instantiation of one or more of the five base classes: Molecule, Basis, QCMethod, Property, and Environment, which are used in the compute and data processing stages. The classes have been designed to be library-agnostic so that they can be used and reused for all libraries. However, not all classes are necessarily needed. It depends on the specific type of calculation that is performed. For example, the first four (or all five if an environment model is used) are needed to run, e.g., a TPA calculation employing Dalton or LSDalton, whereas other functionality may only need some of them (e.g., LSDalton also exposes more fine-grained functionality as illustrated in Section IV A).
The Molecule class contains information about the molecular structure, which can be a single atom, a molecule, a fragment of a molecule, or a set of molecules. It requires, as a minimum, that atomic elements and coordinates are defined. Reasonable defaults are used for other attributes, such as the total charge, spin multiplicity, atomic masses, atomic labels, and, if enabled, information related to point-group symmetry. The atomic elements, coordinates, and, optionally, labels, can be given either as a file, e.g., in XYZ format (optionally with atom labels in the fifth column), or as a string. The DP platform also provides a function that can read the standard molecule file format used by Dalton and LSDalton and return instances of the Molecule and Basis classes. The atomic labels are used in the Basis class as described below but also to specify, e.g., ghost atoms and point charges.
The Basis class holds all basis-set information, which includes the main basis set and, optionally, auxiliary basis sets used for RI and ADMM approximations. The basis set can be given either as a string, in which case the specified basis set is used for all atoms, or as a dictionary using atom labels as the keys and basis sets as the corresponding values. The basis sets are obtained from the basis set exchange Python package. 129 The QCMethod class specifies the method, for example, HF, DFT, MCSCF, CC, etc., together with any associated attributes, such as XC functional or definition of active space.
Approximations used to compute Coulomb and exchange contributions are also given here (and require also that the corresponding basis sets are defined in a Basis instance). In addition, this class is used to specify additional settings, such as convergence thresholds, maximum number of iterations, and so on.
The Property class is used to define which properties to compute. This includes singlepoint energy, geometry optimization, excitation energies, and so on, as well as additional specifications related to the property, which could be, e.g., the algorithm to use in the geometry optimization or the number of states to include in the calculation of excitation energies. Currently, the DP platform supports only a limited set of properties out of the great number that are available in Dalton and LSDalton, but it will be continuously extended. tionality it provides and how the API of the library itself is defined. However, the interface API functions that are exposed to the user must conform to the standards laid out by the DP platform to ensure that there is uniform access to all libraries as well as interoperability between them. This means, for example, that libraries with similar functionality must also provide API functions with identical names and signatures. Moreover, apart from the base classes outlined above, the types and data structures must be either Python built-ins, e.g., integers, floats, lists, dictionaries, etc., or NumPy arrays.
The default ordering on the DP platform is the Dalton ordering: atoms are ordered according to the user input, and for each atom it is angular-momentum components first, and contracted functions second. Other AO orderings may be used on the platform, but the interface API must provide transformation functions to and from the default Dalton ordering. Transformations can then be made on the platform in the cheapest way possible, e.g., on the MO coefficients rather than on the 4-center integrals for MP2, and so forth.
Numerical compute-intensive tasks are performed at the compute stage. This can include anything from the calculation of integral components, all the way to a complete calculation of a molecular property, as well as more complex workflow protocols. In the first development release of the DP platform, users will be able to directly calculate one-and two-electron integrals using LSDalton through CFFI-based Python bindings and have seamless access to a selection of wave functions and properties computed by Dalton and LSDalton using the file I/O communication mechanism. To deal with complex molecular systems that are too large for a full quantum-chemical treatment, we provide an interface to PyFraME that enables workflows involving fragmentation and quantum-classical embedding.
The results obtained at the compute stage can be used directly at the data processing stage, which includes, e.g., data extraction, manipulation, analysis, and visualization. The feature set that will be available in the first development release includes the ability to perform vibrational analysis to obtain vibrational frequencies and normal modes, natural orbital occupation analysis to assist in the selection of active spaces, and partitioning of properties into atomic and interatomic contributions using the LoProp approach. In addition, the DP platform has capabilities for plotting spectra and visualizing orbitals, electrostatic potentials, and electronic densities.
The DP platform can automatically detect and manage available hardware resources.
Manual specification is also possible and allows fine-grained control. In auto-detection mode, the DP platform will first check for resources reserved through a standard HPC job scheduler and, if no scheduler is found, fall back to use resources on the local computer. Resource usage is optimized according to the capabilities of the used libraries, exploiting OpenMP, MPI, or hybrid OpenMP/MPI parallelism when available. Moreover, a task farming functionality is provided for use cases where several separate calculations are to be performed.
We conclude this section with a brief walk-through example that illustrates how the DP platform can be used to streamline a workflow going from initial molecular structures to TPA spectra, employing both Dalton and LSDalton. In Section IV, additional illustrations are presented to demonstrate the fine-grained access to integrals and showcase some of the new features available on the DP platform.
Assuming that the DP platform has been imported as import daltonproject as dp and a list of XYZ filenames is available in xyz_filenames, we start by creating a list of Molecule Then we create the Basis, QCMethod, and Property instances that will be used for a geometry optimization of the molecules using LSDalton. The order in which the classes are instantiated is unimportant. We here start with QCMethod b3lyp = dp.QCMethod(qc_method='DFT', xc_functional='B3LYP', coulomb='DF', exchange='ADMM') specifying that we want to use the B3LYP XC functional together with density-fitting (DF) and ADMM to accelerate the calculation of Coulomb and exchange contributions, respectively. Both acceleration techniques require an auxiliary basis set which is specified when instantiating the Basis class small_basis = dp.Basis(basis='pcseg-1', ri='def2-universal-JKFIT', admm='admm-1') large_basis = dp.Basis(basis='pcseg-2', ri='def2-universal-JKFIT', admm='admm-2') Finally, we create a Property instance geo_opt = dp.Property(geometry_optimization=True) using the default optimization algorithm. With these instances we can proceed to the compute stage. For Dalton and LSDalton, when used through the file I/O mechanism, we use the compute() function. It returns an OutputParser instance that contains methods for transparently fetching relevant results from the output files as shown further below. The compute() function creates a unique filename based on the specific input which is stored in the OutputParser instance. The filename can also be specified through an optional argument in which case it is up to the user to ensure its uniqueness.
We can iterate through the list of molecules, one by one, and perform both the pre-

20
This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. LSDalton can exploit the available CPU resources using a hybrid OpenMP/MPI scheme.

PLEASE CITE THIS ARTICLE AS
The optimal use for LSDalton is typically one MPI process per CPU and one OpenMP thread per CPU core, up to a maximum of 20 threads. By default, additional MPI processes will be used if the number of cores exceeds the maximum. For example, on a computer with two CPUs per node and 12 cores per CPU, LSDalton will use two MPI processes per node and 12 OpenMP threads per MPI process. It is here assumed that the molecules only contain hydrogens, carbons, and oxygens, and that the XYZ files do not contain atomic labels in the fifth column, in which case the atomic labels default to the element symbols. Finally, we create a new Property instance for the TPA cross-section tpa = dp.Property(two_photon_absorption=True) tpa.two_photon_absorption(states=4) where we have additionally specified that we want to include four states. We then proceed with the calculation of the TPA cross-sections using Dalton By default, Dalton will here adopt a purely MPI-parallel execution, corresponding to one MPI process per CPU core. Alternatively, the task farming functionality can be used by supplying the list of molecules to the compute() function tpa_results = dp.dalton.compute (molecule=molecules, basis=mixed_basis, qc_method=camb3lyp, properties=tpa) In this scenario, a separate calculation will run for each molecule, concurrently dividing the available CPU resources amongst them.
Finally, we can plot the TPA spectra using the spectrum module of the DP platform ax = dp.spectrum.plot_two_photon_spectrum (tpa_results) where ax is an instance of Matplotlib 130 Axes class. Further customization can then be performed to produce publication-ready figures.

IV. DP PLATFORM ILLUSTRATIONS
In this section, we provide six use cases of the DP platform with the intent to demonstrate novel platform functionalities in terms of data exposure, processing, and visualization as well as to illustrate some of the added features of the platform libraries. The Python scripts used for these DP platform illustrations, along with the corresponding input/output files, are deposited at https://doi.org/10.5281/zenodo.3710462.

A. NumPy-exposure of one-and two-electron integrals
We here demonstrate how to access primitive and contracted integrals on the DP platform. This functionality is made available by LSDalton through CFFI-based Python bindings and enables access to a variety of one-and two-electron integrals including Coulomb and exchange integral matrices. We will soon add exposure of geometrically differentiated integrals (in principle to arbitrary order), integrals for embedding techniques involving interaction with classical charges and higher-order multipoles, Gaussian-geminal type integrals needed for F12-type theories, attenuated two-electron integrals, multipole-moment integrals, and magnetically-differentiated London integrals (to first order for the two-electron integrals).

22
This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

23
This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.
We note that this example is only meant to illustrate the ease with which algorithms can be implemented on the platform, and that there are better ways to construct the density-fitted  The calculation of UV/Vis spectra with the MC-srDFT method is here illustrated with the retinylidene Schiff base chromophore as originally addressed in previous works. 31,70 Spectra for this system are shown in Fig. 3 together with the π-orbitals that constitute the active space, denoted by π 1 -π 6 . We first note that range-separated calculations introduce a rangeseparation parameter µ that affects the results as long as one remains short of the full-CI

24
This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5144298
limit. From a practical point of view, the optimal µ-value is that which delivers reliable results with the least amount of computational efforts. Benchmark studies have shown that a value around µ = 0.4 a −1 0 is close to optimal for ground and excited states. 26,133,134 We have adopted this value for the present illustration. The following recipe was employed in the spectrum calculations: First, the natural orbitals were calculated from the MP2-srPBE ground-state wave function and the occupation numbers were employed to select a suitable active space. The DP platform provides seamless processing of the MP2-srPBE results and, based on a user-defined selection criterion, an automatic selection of strongly and weakly occupied orbitals is made. Our adopted CAS (6,6) active space came as the result of including orbitals with occupation numbers below 1.99 (occupied space) and above 0.01 (virtual space). Second, the CAS(6,6)-srPBE calculation is performed. Third, the DP platform processes data by extracting excitation energies and oscillator strengths, and generates the requested absorption spectrum. The following definition of the frequency-dependent molar absorption coefficient ε(ω) was used where N A is Avogadro's constant, e is the elementary charge, m e is the electron mass, c is the speed of light, 0 is the vacuum permittivity, n is the refractive index (here set to unity), f i is the calculated oscillator strength of the i th transition, and ω i is the corresponding transition angular frequency. Eq. (3) also introduces the line-broadening function g with the phenomenological parameter γ i corresponding to the half-width at half-maximum (HWHM).
The DP platform offers spectral broadenings based on normalized Gaussian or Lorentzian line profiles.
With outset in Fig. 3, we briefly discuss some key features of the MC-srDFT method. We note that the CAS-srPBE spectrum shows two distinct peaks in good agreeement with the experimentally observed S 1 and S 2 bands 135,136 -transition energies and oscillator strengths are provided in Table I. It is here essential to employ a multi-configurational wave function as can be seen by the comparison to the corresponding single-determinant (HF-srPBE) spectrum, which essentially corresponds to the range-separated LC-PBE model. Without account of static correlation, peak positions are strongly blue-shifted and the intensity of the second band is severely underestimated. It is noted that occupation numbers from a regular MP2 calculation suggest significantly larger active spaces. 137 However, the use of CAS(6,6)-srPBE provides the same accuracy as literature CASPT2 results based on a much larger CAS (12,12) active space. 138 This demonstrated opportunity to employ relatively small active spaces with MC-srDFT is not limited to the retinylidene Schiff base chromophore but has also been shown in other contexts, e.g., to describe the mechanism of [NiFe]-hydrogenase. 139

26
This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5144298
The reason behind the failure of single-reference DFT approaches to properly describe the electronic transition underlying the S 2 band in the retinylidene Schiff base chromophore becomes apparent from an analysis of the MC-srDFT response vectors. For S 1 , the dominant element of the response vector refers to the configuration associated with a single-electron excitation from π 3 to π 4 (CI coefficient of 0.78 in Table I). For S 2 , on the other hand, there are three almost equally important configurations appearing in the response vector and one of which is associated with a double-electron excitation from π 3 to π 4 (CI coefficient of 0.42 in Table I) that is unaccounted for in standard time-dependent DFT.

C. Modeling complex systems through fragment-based quantum-classical approaches
We here provide an illustration that demonstrates not only the ability of the DP libraries to perform spectrum simulations of complex systems but also the ability of and potential for the hosting DP platform to manage the workflow of complex computational protocols associated with environment fragmentation and configuration-space sampling. It is an indisputable fact that first-principles methods in quantum chemistry come with a computational cost that hampers applications to relevant chemical and biochemical systems such as solutions and protein-embedded chromophores. Methods that allow for low-cost approximate yet accurate modeling of environments are therefore scientifically enabling, and one such approach is the PE model briefly described in Section II B 3. In the PE scheme, the environment is represented by atom-centered multipoles (typically up to and including quadrupoles) and atom-centered anisotropic dipole-dipole polarizabilities that allow for a mutual polarization between the quantum part and the classical environment.
As further described in Section II B 3, the adopted multipole expansion in the standard Both the PE and PDE models are applicable to large and complex (bio-)molecular systems. Since the parameters describing the environment are derived based on a fully ab initio description, the setup of the spectrum calculation involves a large number of preliminary calculations on the separate fragments of the environment. Figure 4 illustrates this situation in terms of TPA spectrum calculations of the Nile red chromophore embedded in the β-lactoglobulin protein. Including the solvent, this system contains 32,582 classical sites.
The entire protein and water molecules within 15 Å of the Nile red chromophore were treated with either PE or PDE while water molecules beyond this distance were described using the TIP3P water model (in total 10,000 water molecules).
Configuration sampling of the environment was performed by averaging over 150 snapshots sampled from a quantum mechanics/molecular mechanics (QM/MM) MD simulation. 140 The calculation of embedding parameters for each snapshot took circa 1 and 4 hours for PE and PDE, respectively, on 20 nodes (dual socket E5-2680v3, 12 cores). The TPA spectrum calculations (five lowest singlet states) included EEF effects and took circa 4.5 and 5.0 hours per snapshot for PE and PDE, respectively, on 8 nodes.
Assuming a monochromatic laser source, the TPA cross-section (in units of GM) is computed as 141 where α is the fine-structure constant, a 0 is the Bohr radius, and g is here chosen to be a Lorentzian with a HWHM of γ = 0.1 eV (see further Eq. (3)). The two-photon (TP) transition strength of the i th transition is computed assuming linearly polarized light as where S ab are the associated TP transition matrix elements.
The gas phase TPA spectrum of Nile red shows three distinct peaks at 1.5, 1.9, and 2.2 eV. Embedded in the protein, the lowest band becomes red-shifted by 0.1 eV and it also becomes both broader and lower in intensity. At higher energies, the PDE model preserves the structure from the gas phase partially, with an intense peak at 2.1 eV, and a shoulder at 1.7 eV. In contrast, the PE model predicts a large peak at 2.1 eV with a broad shoulder towards lower energies.
D. Open-ended response theory for electromagnetic and geometric perturbations: Infrared and Raman spectroscopy with the OpenRSP and

SpectroscPy modules
Spectroscopic techniques involving molecular vibrations are useful for characterizing molecular systems, and with the DP platform, an option for the calculation of vibrational spectra is available through the combined use of LSDalton, OpenRSP, 120 and SpectroscPy. 142 We here illustrate this functionality through a calculation of IR and Raman spectra of benzene, including a brief outline of the key aspects of the LSDalton/OpenRSP/SpectroscPy combination and its use on the platform.
From the DP platform, energy-derivative tensors generated by LSDalton/OpenRSP are passed on to SpectroscPy that processes them to generate spectroscopic properties.
OpenRSP manages the calculation of response properties by an open-ended quasienergybased formulation of response theory 143 that employs a recursive formalism. 144 OpenRSP has an API through which it can be connected to different host programs. OpenRSP is currently called through LSDalton so as to provide the necessary integral derivatives (briefly outlined in Section II C 1), XC contributions through modules XCFun 95 and XCint, 121 and response equation solver capability. 145 The use of QcMatrix 128 allows OpenRSP to be agnostic to the details of the matrix operations, i.e. independent of their specific implementation when passing from/to the host program and when carrying out operations inside the OpenRSP core functionality. SpectroscPy is used to produce spectroscopic properties involving molecular vibrations, by performing vibrational analysis generating vibrational frequencies and absorption properties. Presently, SpectroscPy offers functionality for IR and Raman spectroscopy in the harmonic approximation, and can also combine data from calculations on a series of molecular configurations, which is useful, for example, in applications dealing with flexible molecular systems that require configuration-space sampling. Future extensions will include anharmonic corrections, hyper-Raman spectroscopy, and other spectroscopic processes.
The IR molar decadic absorption coefficient for normal mode i can, in the double harmonic approximation, be expressed as 19 where µ is the molecular dipole moment, Q i is the normal mode coordinate i, and g is here Similarly, the harmonic differential Raman scattering cross-section σ i is given by 146 where h is Planck's constant,ν 0 the wavenumber of the incident light, k the Boltzmann constant, and T the temperature. We have here introduced where α is the molecular dipole-dipole polarizability.
The spectra shown in Fig. 5 were generated by SpectroscPy using the molecular Hessian and the first-order derivatives of the molecular dipole moment and polarizability calculated This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS
spin-restricted formalism ensures that the ground-state electron density is free from spincontamination and, as such, it provides a better starting point for molecular property calculations compared to the more widely used unrestricted Kohn-Sham (UKS) approach. The advantages of ROKS over UKS are most apparent in calculations of electron paramagnetic resonance (EPR) spin Hamiltonian parameters, which are explicitly dependent on an effective expectation value of the electronic spin operator, and benchmark studies on organic radicals show that the ROKS approach is able to better predict electronic g-tensors and hyperfine coupling constants. [151][152][153][154] The main strength of the ROKS approach is the ability to produce spin-contamination free electron density for the high-spin ground state. This is achieved by imposing spin-symmetry restrictions during the SCF optimization. 149 The side effect of these restrictions is that the hereby obtained KS orbitals have a non-vanishing gradient with respect to orbital rotations of triplet spin-symmetry, 152,155 and this side effect manifests itself in the computation of electronic spin-dependent properties, like hyperfine coupling constants. For example, the expectation value of the one-electron spin-dependent operatorÂ in the ROKS approach is given by the direct spin-density and spin-polarization contributions 152 Â = Tr (AD spin ) + Tr (AD pol ) ; where the first contribution is computed in the same way as in the UKS approach by contracting the electron spin-density D spin with operator matrix A in AO basis, and the second contribution is computed similarly to the first contribution by replacing D spin with the spinpolarization density D pol .
The spin-polarization density, D pol , is determined by solving restricted-unrestricted response equations that account for the relaxation of KS orbitals in the presence of the spindependent perturbation. The relative importance of the spin-polarization density contribution varies greatly between different molecular properties: from being prominent for isotropic hyperfine coupling constants 152,155 to being negligible for electronic g-tensors. 154 The use of the spin and spin-polarization densities as obtained from the restrictedunrestricted response formalism is not limited to the calculation of spin-dependent molecular properties but it can also be employed for topological analyses of open-shell system responses to spin-dependent perturbations. To illustrate this point, we depict isodensity surfaces of the spin and spin-polarization in Fig. 6 for the 2 B 1g ground state of copper acetylacetonate This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. complex, Cu(AcAc) 2 . The spin-density contribution for Cu(II) complexes is expected to primarily stem from the 3d xy -orbital of the copper atom, 154 but it is here seen to also acquire significant contributions from the coordinating oxygen atoms. The spin-polarization density is also delocalized across copper and oxygen atoms and based on the localization of these two densities, one can predict the behavior of spin-dependent properties and qualitatively estimate the importance of the spin-density and the spin-polarization contributions. We emphasize that the presented disentangled visualization of spin-density and spin-polarization contributions to molecular properties is uniquely accessible in the restricted-unrestricted response formalism and not available in the UKS approach.

F. Coupled cluster methods for inner-shell spectroscopy
As an illustration of the use of the DP platform for the simulation of X-ray spectroscopies, we present in Fig. 7  The three lowest excitations at the near carbon edge correspond to the first two bands in the experimental spectrum, These excitations are assigned to 1s → π * transitions from the three carbon atoms by means of a natural transition orbital (NTO) analysis. The transition associated with the carbonyl group is chemically blue-shifted by some 1.5 eV from the other two, see Fig. 7. There is an overall excellent spectrum agreement between theory and experiment and, although simple, this example illustrates well the value of spectrum simulations for the characterization and interpretation of experiments. The most prominent added value of the DP platform for ground-state X-ray spectroscopy simulations using Dalton comes at this stage of analysis and visualization.
In the example above, the CVS approximation is employed only during the determination of the core-excited states. It is implemented as a projector that during the iterative solution of the CC eigenvalue problem, 56,57 AR f = ω f R f , only retains elements R f aI , R f aI,bj , R f ai,bJ ,

34
This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. and R f aI,bJ (CCSD case), where I and J refer to core orbitals and a and b refer to virtual orbitals. A schematic representation of the CVS approximation is shown in Fig. 8, where only the green matrix sub-blocks are effectively kept after applying the projector during the iterative procedure of the response solver. Core and valence excitations become decoupled and X-ray spectra are obtained at basically the same computational cost as UV/Vis spectra with use of identical bottom-up iterative techniques. The ability of a specific CC method to reproduce specific core-excitation spectral signatures depends on the amount of relaxation as well as of single, double, or higher excitation character of the transition. CCSD generally yields core spectra in rather good (semi-quantitative) agreement with experiment, with rigid blue-shifts ranging in between 0.5 and 3 eV, depending on the K-edge considered and the basis set adopted.

PLEASE CITE THIS ARTICLE AS
In addition to facilitating data analysis and visualization, the DP platform greatly simplifies the more complicated setup involved with simulations of transient (excited-state) X-ray absorption and emission spectra. Due to design legacies in Dalton, such simulations without the use of the DP platform will require a substantial amount of manual file handling as the necessary valence and core excitation vectors must be obtained from separate code executions and then later recombined. On the platform, all of these steps are handled in an automated way.

V. OUTLOOK
We have given a presentation of the Dalton Project that marks a paradigm change in the software engineering practices for the Dalton ecosystem. At the heart of the DP, we find a hardware-aware platform written in Python with support from NumPy, SciPy, and MPI4Py which provides a modern user interface and defines a communication standard for seamless library access and interoperability. At present, the DP platform supports three libraries namely Dalton, LSDalton, and PyFraME but further modular extensions are underway.
Libraries are written in any of the programming languages Python, Fortran, C, or C++ and they communicate with the platform by means of file I/O, Python bindings though CFFI or pybind11, or direct Python module import. Modular programming with enforced unit testing will become the standard for newly started developments and it is to be anticipated that the existing few monolithic codes will gradually be phased out and replaced by a larger number of libraries with more specific tasks. The DP platform is open source and distributed under the GNU General Public License version 3.0 or later (GPLv3+).
The DP platform is developed for execution on personal computers as well as more powerful servers in HPC environments and, although without guarantees, the user should expect it to run under Windows, MacOS, and Linux operating systems equipped with Python 3.6 (or later) installations, see the DP website https://daltonproject.org for instructions.
Installation of the DP libraries is a separate issue and the platform will gracefully signal to the user when it encounters called for but missing libraries. It is of course fully possible to use the DP platform alone on a personal computer to analyze the results available in output files produced on a remote HPC system. As our model to reach a sustainable software ecosystem adopts the notion of separate and quite independent release and distribution policies for the libraries, the main burden of work and the most dependency issues in the installation process are expected to be found in the phase of library installation. Driven by the stimulus of becoming recognized and used, it is anticipated and expected that newly started library developments will care to offer a smooth installation process on a widespread selection of platforms and operating systems.
For developers, the DP platform can lead to rapid prototyping of novel scientific ideas as illustrated in Section IV A with an examination of the RI approximation. But this example also points out something else of great importance namely the educational aspect of the 36 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. PLEASE CITE THIS ARTICLE AS DOI:10.1063/1.5144298 DP platform. Our experience tells us that the process of implementing methods to solve fundamental equations is supremely efficient to reach a deeper understanding of the topic at hand, but only few students are granted this opportunity as core program modules of scientific software are written a long time ago and often made obscure by code optimization.
What is here illustrated is the access to the needed building blocks to explore quantum chemistry in very much the same manner that we can use the Python NumPy package to explore linear algebra. At the early stage of a PhD education, we believe that this can be of high value and help overcome some of the initial hurdles faced during a career in quantum chemical theory and program development.