Jellyfish: A modular code for wave function‐based electron dynamics simulations and visualizations on traditional and quantum compute architectures

Ultrafast electron dynamics have made rapid progress in the last few years. With Jellyfish, we now introduce a program suite that enables to perform the entire workflow of an electron‐dynamics simulation. The modular program architecture offers a flexible combination of different propagators, Hamiltonians, basis sets, and more. Jellyfish can be operated by a graphical user interface, which makes it easy to get started for nonspecialist users and gives experienced users a clear overview of the entire functionality. The temporal evolution of a wave function can currently be executed in the time‐dependent configuration interaction method (TDCI) formalism, however, a plugin system facilitates the expansion to other methods and tools without requiring in‐depth knowledge of the program. Currently developed plugins allow to include results from conventional electronic structure calculations as well as the usage and extension of quantum‐compute algorithms for electron dynamics. We present the capabilities of Jellyfish on three examples to showcase the simulation and analysis of light‐driven correlated electron dynamics. The implemented visualization of various densities enables an efficient and detailed analysis for the long‐standing quest of the electron–hole pair formation.

into the attosecond timescale are possible. 1Such experiments help to understand the light-matter interaction in great detail for a variety of materials ranging from small molecules to nanostructures.Ultrafast electron dynamics plays a central role in the intricate ionization processes with strong fields of molecules, that is, high-harmonics generation (HHG), 2,3 charge migration processes in organic molecules, [4][5][6] charge transfer in heterogeneous photocatalytic processes, 7 or energy transfer in semiconductor nanostructures. 8,9In many studies of photo-induced processes in biological systems, few electronic states are considered with the focus on nuclear dynamics computations and the involvement of conical intersections, which play a significant role in the UV damage-and protection-mechanisms of DNA, 10,11 for example.The photochemical properties, however, are inherently the result of a state manifold, detailed understanding of radiation-induced electronic mechanisms, and the determination of excitation time scales and lifetimes of the states are of central interest.
Also in inorganic materials, the investigation of electronic processes in the time domain experiences a rising interest.One famous example is semiconductor quantum dots (QDs) that are increasingly gaining in importance for the usage as qubits in quantum computers, [12][13][14] for example.For the understanding of the capabilities of QDs as qubits, dynamical manipulation processes are investigated. 15Moreover, undesired side effects are studied in order to estimate and suppress their influence. 16,17These can be, among others, the Auger decay within a QD [18][19][20] or the inter-Coulombic decay process between two or more QDs. 21,22Both processes lead to the ionization of a QD and may decohere desired qubit states.
Even if many dynamical properties of a system like ionization rates or absorption spectra can be predicted from time-independent properties like transition dipole moments or excitation energies, the faster the considered processes happen and the more electronically excited states are involved, the more it is necessary to simulate the processes explicitly and to solve the time-dependent Schrödinger equation for the considered time period.
The methods used for this purpose are often extensions of the corresponding time-independent theories.There are single-reference methods, like time-dependent Hartree-Fock (TDHF) 23,24 or time-dependent configuration interaction (TDCI) [25][26][27] to propagate a wave function in time explicitly.Also, time-dependent methodologies that propagate a reduced density matrix ρ in time to include effects of dissipation exist. 28,29Multi-reference time-dependent methods like multiconfiguration time-dependent Hartree-Fock (MCTDHF) 30,31 are needed for situations with high static correlation like a photo-induced dissociation, for example.Even for the static density-matrix renormalization group (DMRG) methodology 32,33 there is a possibility for real-time evolutions by combining it with time-evolving block decimation (TEBD). 34or electronic structure calculations, various programs like ORCA, 35 Psi4, 36 Q-Chem, 37 PySCF, 38 and many more are available for academic and commercial users, which were developed over decades by large user and developer communities.They cover all steps from reading in and optimizing molecular structure to determining properties like nuclear magnetic resonance shielding parameters, dipole moments, vibrational, and optical spectra at highest efficiency.This is not the case for most codes in the field of electron dynamics methods due to the younger age of the research field and the relatively smaller research community.Hence, available code is currently limited to executing the essential dynamic steps.This means that at the current state, to perform and illustrate such dynamics calculations, several programs are needed and intermediate results have to be handled through various input and output files.A possible procedure for TDCI calculations, as done for example by Weber et al., 39 consists of four steps.After geometry optimization, the molecule and basis set are read into ORCA to calculate the ground and all excited states using CI.The ORCA output file is then used in ORBKIT [40][41][42] to determine the full dipole matrix.Based on the ORCA output file and the dipole matrix, a ρ-TDCI calculation is performed using the GLOCT program. 29Finally, the time-dependent wave function is analyzed and properties of interest are extracted.As can be seen from the example, many programs or tools are often designed for a single purpose-a single propagation method or a certain use case-and little attention was given to the possibility of successively extending them.
In this work, we introduce a platform with the program Jellyfish, which carries out electron dynamics simulations starting from the electronic structure calculation and ending with the visualization of densities in motion in a single program.The main objective of the Jellyfish package is to provide a single program suite that is most easily expandable without becoming complicated for developers and users.These targets are obtained by a module and a plugin system.Individual partial calculations such as the calculation of dipole moments, the reading-in of a molecular geometry, for example, are developed as separate modules.These modules use the results of other modules and additional parameters for the successive calculations and then provide output for the following modules.In this way, flexible computational workflows can be assembled by modules for a large variety of application purposes.Uncomplicated creation and modification of a workflow is possible through the graphical interface (GUI) that employs a NodeEditor to logically assemble and connect the modules.The incorporation of such a NodeEditor was inspired by other fields like 3D design, where the use of such a system has become a common interface.In combination with the use of a GUI instead of solely textbased input and output files allows an intuitive way of using the program.This makes it especially easy for new users to get started with the dynamics simulations.
Furthermore, a plugin system is underlying the code that allows the development and compilation modules independent from the main Jellyfish code.Several modules are bundled into a plugin and dynamically loaded into Jellyfish through an interface.New or modified methods and modules can thus be easily exchanged and published as plugins.To simplify the exchange of modules and the further development of Jellyfish, the code at its present state as well as the existing plugins, are published under an open source license (GPL-3.0)on github: https://github.com/FabianLangkabel/Jellyfish.The plugins developed so far have focused on the TDCI method, the visualization of time-independent and time-dependent wave functions, densities, and orbitals, as well as fundamental quantum computing algorithms.With the functions available so far, Jellyfish can thus be used immediately for a number of dynamics problems without further development.
The article is divided into four sections.In Section 2, the exact requirements for Jellyfish are formulated and it is explained in how far these requirements are met by the design of the program.Here, the already mentioned module and plugin systems, as well as the graphical user interface and the file format, are discussed.The Section 3 discusses which functions in the form of modules and plugins are already available with the release of the program.Section 4 deals with three examples, which have been computed and evaluated with Jellyfish and demonstrate the functionality and flexibility of the program as well as key functions of the already implemented methods.In Section 4.3, the main aspects and possibilities of the presented program are summarized.

| DESIGN AND IMPLEMENTATION
The main objective of Jellyfish is to provide a platform that performs all steps necessary for a TDCI calculation and its visualization in a single program suite.To do so, a chain of calculations and their postprocessing can be set up in a consistent manner so that the users do not have to handle several programs and formats.That way, the margin for errors is reduced and, for example, conflicts due to updates for the data transfer from one program to the next can be avoided.In addition, the following demands were made for the development of Jellyfish: • Flexible: Jellyfish aims to provide a general interface in which various Hamiltonian, basis sets, wave function analysis tools, and other partial calculations can be arbitrarily combined as desired.For example, it should be possible to perform dynamics not only in molecules, but also in a general electron-binding model potential.Also, a variety of propagators and time-dependent methods exist besides TDCI simulations.The time-dependent Hamiltonian also varies from application to application and may include, for example, interaction with external laser fields to compute laserdriven dynamics or complex absorption potentials to simulate ionization processes.

• Extensible:
The open-source Jellyfish code should enable users to develop their own modules or routines.If possible, the development of extensions should not require in-depth knowledge of the program architecture and no or as little as possible modification of the existing program code.• User-and developer-friendly: Jellyfish should be user-friendly and self-explanatory.The goal is to keep the barrier to using the program and the development of new code as low as possible.• Interactive: Often simulations require the evaluation of prior calculations and therefore cannot be automated.For example, for light-driven simulations, the knowledge of target state properties is needed to determine the parameters of an external field.Jellyfish should enable such an analysis with subsequent adjustment and evaluation of the calculation on the fly, without the need to close the program and save and read in data in the form of various input and output files.• Efficient: Like almost all programs in the field of quantum theory, Jellyfish should run as efficiently as possible to allow for fast processing and the largest possible systems.Not only should it be written as an efficient code, but it should also enable to run calculations on headless servers and in combination with common parallelization options like MP, MPI, or CUDA.
The C++ programming language is chosen for these purposes as it is inherently object oriented and widely used.C++ allows the compilation of very efficient code, for which a large number of libraries are available.Especially the application framework Qt which was used for the development of the here presented user interface was crucial.

| NodeEditor/visual dataflow programming
In order to meet the above requirements in an optimal way, Jellyfish orients on dataflow programming.In dataflow programming, a program is modeled as a directed graph in which data runs through various operations and is modified by them.This approach is realized in this work by a graphical user interface with a NodeEditor.Such approaches were already successfully used in programs from other fields, for example, in 3D and game design.Each partial calculation or execution-may it be the reading in of a molecular geometry or carrying out an electronic structure calculation-is not only developed as a module, but also represented graphically in a network in the NodeEditor by an individual node in form of a box with input and output ports, as displayed in Figure 2. The implementation of such a partial calculation and everything associated with it is called a module, and instances of a module are called nodes.In a project, there can be several nodes instantiated from a module.For example, two nodes for generating basis sets.The nodes have different identifiers to distinguish them.The network that is built from nodes is also referred to as a project in the following.
The input ports on the left side of a node box represent data types, for example, a matrix, and the output ports on the right of a node can be connected arbitrarily with the input ports of other nodes, as long as they have the same data type.The data are then propagated from one node to another and shared between the connected nodes.The node that creates the data (output node) always holds the data, while the nodes that use the data as input receive a reference (pointer) to that data.Such an approach allows the flexible construction of individual calculations for the most flexible purposes and can be extended or connected arbitrarily to other or new modules without having to modify existing code.Furthermore, the graphical approach avoids complicated input formats so that Jellyfish can be handled in a much more self-explanatory way than is possible by pure command line interfaces.Example networks for various simulations are published along with the code.Since a network of nodes can be extended and modified arbitrarily, even if other nodes are already calculated, calculations can be performed and modified interactively step-by-step.Furthermore, a new node may be added to the network from any type of source code in form of a plugin, see next subsection.
The underlying definitions of modules, nodes, connection, and project files, as well as the related functionalities, have been developed as JellyfishCore.JellyfishCore is compiled as a static library and forms the foundation for the Jellyfish architecture.The schematic structure of Jellyfish and the dependencies of the individual components is shown in Figure 1.

| Plugins
To separate the code of the modules from the user interface and to simplify the development of new modules, several modules are combined into one plugin.For example, different analysis modules were combined into the analysis plugin.The compilation of the plugins is independent of the user interface.The plugins are based on the JellyfishCore and contain the NodeRegistry, which creates, destroys, and stores instances of previously defined modules.Additionally, they have an interface that allows communication with the main program.Plugins are compiled as shared libraries, then placed in a plugin folder of the Jellyfish installation, dynamically loaded and included.As a blueprint, the C++ library micoplugins is used and heavily modified in-house for our specific application purposes.This plugin structure allows modules or entire plugins to be easily added without having to modify or recompile the other program parts.Plugins can thus be developed separately and easily exchanged.

| JellyfishGUI
JellyfishGUI is the primary input program with its graphical interface and a screenshot of the main window is shown in Figure 2. On the left side, it offers a list including a search function, which displays all available modules and offers the possibility to add them as a node to the network via drag and drop to the central part, the NodeEditor.In this space, see Figure 2, the connecting and disconnecting of NodePorts is done.To include a NodeEditor in the Qt-based user interface the library QtNodes 43 is adapted to the JellyfishCore and included in the code.
By selecting a node in the NodeEditor, node-dependent details and setting options are displayed on the right-hand side.In the example shown in Figure 2, the propagation module was selected, which is why propagation parameters such as the number of time steps and their size can be set accordingly.Also, further functions or whole windows for individual nodes can be integrated.In the case of the propagation module, for example, a whole window is available behind the "Analyze Propagation" button in the figure, which is described in more detail in Section 4.1.Other examples for such additional functions are the export of a laser pulse as ASCII file or the integration of the subprogram for the analysis of time-dependent electron densities as it was done in the Visualization plugin.Furthermore, the main window, shown in Figure 2, offers the view of a project log, as well as the possibility to work on several projects at the same time by tabs and, of course, the loading and saving of project files.

| JellyfishCMD
For the usage on headless servers in an HPC environment, JellyfishCMD was developed.It is a command-line interface allowing the calculation on individual nodes or entire node networks and, therefore, can be managed through scripts.A typical procedure for expensive computations in Jellyfish is the creation of a project in JellyfishGUI, the execution of the calculation on an HPC with JellyfishCMD, and subsequent evaluation with JellyfishGUI.

| File format
All data of a project are stored in a single Jellyfish file (file name extension .jlf).It is a ZIP archive created and read by employing the zlib and libzip libraries.Besides the computed data objects, for example, excited state energies, or timedependent populations and alike, this archive also contains module files that store information about the nodes and their connections in a JSON (JavaScript Object Notation, using the JSON for Modern C++ library 44 ) format.In the JSON file, the type of node and its unique GUID (Globally Unique Identifier, generated by the library stduuid 45 ) as well as its status-computed, not computed-are stored.To avoid multiple storages in memory, the data always belongs to the node in which it was calculated, while all connected nodes are assigned pointers to this data.Computed results are stored in files composed of the node GUID and a name of the result, "DensityMatrix" for example.All these data are also written directly into the JLF archive and read from it.The form in which the data are stored can be selected individually in each plugin.Saving the data as a binary file to save storage and to enable faster saving and loading is conceivable.The focus, however, was placed on being able to read the data as easily as possible with external programs.Thus, it is stored accordingly in plain text.The choice of a ZIP archive as the format for the JLF files has on the one hand the advantage that it accounts for less storage space due to the compression that ZIP archives provide.On the other hand, objects are stored in separate, unique files in the archive, which simplifies the external access through scripts and programs.For example, the module file can be accessed and modified by a Python or shell script for further postprocessing.

| Basic plugin
The goal of the Basic Plugin is to provide all functionalities needed to prepare TDCI simulations.For such simulations, various modules have been implemented, starting from reading geometries of molecules or general binding potentials, through a Hartree-Fock section with its evaluation of one-and two-electron integrals and the TDCI calculation, to subsequent propagation of the wave function.
The first module allows the reading of molecular geometries in Cartesian coordinates (XYZ file format), the selection of a basis set and the assignment of a charge.Spherical as well as Cartesian Gaussian basis sets can be selected.The most common Pople and Dunning standard basis sets are readily implemented and more can be downloaded from Basis Set Exchange 46 in the same format as for Gaussian or Psi4.However, the addition of further basis functions with flexible, custom-made parameters is also possible.For example, it is possible to add any arbitrary set of basis functions, for example, an even-tempered basis set as used for the example in Section 4.2.
For the calculation of various one-and two-electron integrals, the library libcint 47 was included in the plugin.Further integrals between basis functions and arbitrary potentials, which are not analytically available, can be computed numerically on a grid and added to the total Hamiltonian.For this purpose, the library muparser 48 was embedded, which allows the parsing of input functions.
To calculate the electronic structure, a restricted Hartree-Fock (RHF) and a configuration interaction module are available.The CI module is based on Slater determinants and can perform calculations at any CI truncation level.Also available are modules for the necessary orbital transformations from atomic orbitals to molecular and spin orbitals, respectively.The implemented electronic structure code is kept as simple as possible and is currently still inefficient, which makes it suitable only for small calculations.As an alternative, however, an interface to the ORCA program package is provided.In a single module, it reads geometry, basis set, Hartree-Fock results, and CI-or TDDFT-states from an ORCA log file and can be used for further processing, for example, for the visualization of orbitals and stationary states or for time-evolution, as demonstrated in Section 4.1.The development of a plugin that integrates other external electronic structure programs directly is planned.The open-source program package Psi4, for example, over the supplied application programming interface (API) for C++, would be suitable as an ideal interface.
As outlined before, subsequent to an electronic structure calculation, correlated electron dynamics simulations by TDCI can be carried out and visualized by corresponding modules.For the time evolution, a module is currently available that propagates the system on the basis of the CI states in time steps using the split operator technique.
For the set-up of a propagation, several options are available.The propagation module offers the choice for the initial state, that is, the ground state and simultaneously any excited state, as well as the possibility to limit the number of CI states used during propagation.In addition to the time-independent Hamiltonian from the previous steady-state CI calculation, other time-independent operators can be added to the operator, both Hermitian or also non-Hermitian ones, where an example for the latter is a complex absorbing potential (CAP).In the case of the CAP, two options are available.A spatial, spherical CAP, which is defined by a radius around the center of the system, 49 or an energetic CAP, which is defined by the energy of virtual molecular orbitals.
Finally, time-dependent simulations are possible through the inclusion of the propagation module and a laser module.The network of a TDCI propagation when starting from an ORCA calculation is displayed in Figure 2. The interaction with a laser field in TDCI is treated in the semi-classical dipole approximation The laser module allows the setting-up of an external electric field F ! t ð Þ and an input-connection for the dipole matrix of the CI states, For the calculation and display of the dipole matrix, another module is implemented, which also uses the libcint 47 library for integration.
The modulation of a laser pulse is also available as an individual module.The pulses have the form where s t ð Þ defines the envelope function and r t, ω ð Þ the driving field function of the pulse.Variations for these functions are readily available.In addition to the amplitude f 0 and frequency ω of the pulse, the pulse duration can be specified in atomic time units or in terms of the number of laser cycles.As polarization, either a linear polarization with an arbitrary direction vector or a left or right rotating polarization can be used.The defined laser pulse can be visualized directly in the GUI before it is used in the propagation.More complicated pulses that cannot be modulated with these options can be imported from an externally generated ASCII file.After a propagation, properties like the populations, the norm of the wave function, the time-dependent dipole moment, or a Fourier transform of the latter to obtain, for example, a HHG spectrum, can be displayed and saved in the propagation module.
In many functions, matrices or tensors are used in the program.For these and most operations from linear algebra, like matrix multiplications, eigenvalue determination, and so on, the library Eigen 50 was used.Eigen has very efficient implementations of the functions, most of them also including parallelization via OpenMP.

| Read external plugin
The aim of Jellyfish is not to provide "another" HF/post-HF program suite for electronic structure calculations but rather to enable their postprocessing with visualizations and/or time-dependent simulations.Even if the Basic Plugins provide these routines for electronic structure calculations, they are, as already mentioned, not yet optimized and thus comparatively slow.
The ORCA program allows such calculation in a more efficient form and a Read External Plugin that reads the log file is included in the Jellyfish release.That way also methods that are not included in the Basic Plugin can be used for Jellyfish processing.For example, it is possible to carry out light-driven dynamics simulations with Jellyfish by importing configuration interaction singles (doubles) (CIS(D)) or linear-response time-dependent density functional theory in Tamm-Dancoff approximation (LR-TDDFT) calculations, as presented in the example in Section 4. Both methodologies improve on the system energetics while retaining the simple form of the CIS wave functions.In CIS(D), the double excitations are included as perturbative corrections 51 to the ground and excited state energies via an energy correction, while the wave functions are still the CIS wave functions.The CIS(D) method has the advantage of working with a CI matrix of the same size as CIS, while at the same time accounting for an improved treatment of electron correlation.In LR-TDDFT, correlation is incorporated through the use of the exchange potential and excited states are obtained by orbital transitions with the ground state configuration as a single-reference.This makes it possible to construct a CIS-like pseudo wave function which then can be used in the TDCI methodology.The Read External Plugin has therefore a module, which makes the import of basis set, molecule geometry, Hartree-Fock or Kohn-Sham orbitals as well as CI vectors and energies from an ORCA log file possible.

| Visualization plugin
A special focus in the program development was given to the visualization plugin.The goal was to develop modules that enable visualization of wave functions and electron densities with a focus on laser-driven electron dynamics.Visualizations are implemented in two modules, one for stationary and one for time-dependent orbitals and densities and properties as listed in Table 1.As a basis for the visualization, the VTK library 52 has been directly linked into the plugin and a window to display orbitals and densities in terms of 3D isosurfaces has been developed, which can be called directly from the JellyfishGUI.
These plot modules offer the choice to also display the molecular structure (stick and ball) as well as a coordinate system.Besides the display of static atomic orbitals χ r 53,54 for both static and time-dependent situations, as summarized in Table 1.
The modules for the visualizations of static and time-dependent properties have a similarly structured window.A screenshot of the developed visualization front end, Plot Orbitals and States, is shown in Figure 3. Here, the timedependent difference density for the nucleobase guanine is shown as an example.The visualization windows show four areas.In the central panel (a) wave function objects, like orbitals or densities are displayed and can be viewed in free rotation controlled by the mouse or keyboard.On the left side of the window (b), the object to be plotted/rendered can be selected (cf., Table 1, depending on the selection further options may become available).Along with a choice of the object type, certain time step or a time interval can be specified.The control of the graphics settings is on the right side of the window (c).Here, the object colors can be chosen along with graphics settings that control the quality, resolution, and speed, as well as options to display additional objects like the molecular structure, a xyz axis, and the display of a time stamp.To enable viewing a video "on-the-fly" in the main display window even for larger polyatomic molecules, a radius to truncate basis functions was introduced.This section also has options to export the in (a) displayed objects/ arrangement as a frame or as a movie.A log file (d) showing updates on the rendering and export progress is given below the main plotting window.
The spacial extend of orbitals and densities is represented by an isosurface composed of the sum of the basis function contributions.Their positive and negative values are displayed in two different colors.To additionally highlight their compactness or diffuseness, an option to also show multiple isosurfaces is available.
T A B L E 1 List of available visualization options for orbitals and electronic state properties.

Property
Formula Stationary Time-dependent

| Quantum computing plugin
Additionally to the above outlined traditional way of simulating TDCI dynamics, Jellyfish also was used to develop algorithms for quantum computing. 55Quantum computing promises for great advancement in electronic structure theory [56][57][58] as well as quantum dynamics. 59Likewise, the field of electron dynamics would benefit from quantum computers and would make exact calculations in much larger molecules finally feasible.The quantum computing plugin for Jellyfish introduces all necessary modules to perform electron dynamics calculations on a quantum computer simulator.For gate-based quantum algorithms, the library QuEST, 60 which provides highly parallelized code, is embedded.Even noisy quantum computers can be simulated with QuEST.Modules available in the quantum computing plugin include the Jordan-Wigner 61 transformation of the time-independent operator and laser fields to translate the integrals of molecular orbitals into weighted Pauli strings for the quantum computer.For propagation, a module for Hamiltonian simulation based on Trotter decomposition 62,63 is available, as well as an extension with the quantum imaginary time evolution (QITE) 64 algorithm to simulate ionization as well.Time-dependent expectation values can be determined with the help of a module for the Hadamard test. 65More information about the complete algorithm developed with it can be found in Langkabel et al. 55 4 | EXAMPLE APPLICATIONS

| Electronic transition in guanine
As a first example to illustrate the capabilities of the Jellyfish package in a many-electron simulation, a photoexcitation of the nucleobase guanine is presented.Guanine (2-amino-6-oxopurine) is an aromatic heterocycle consisting of a 6-membered and a 5-membered ring and has a total of 76 electrons.As a building block of DNA, it has been part of numerous photochemical studies.Highly accurate and computationally demanding multiconfigurational and multireference electronic structure calculations were employed to compute the (minimum energy) reaction paths, for example, to study the efficiency of radiation-free decay mechanisms. 66A comprehensive study on guanine tautomers, for example, discusses spectra and geometries employing CIS and CASSCF. 67There, also electron density difference plots have been employed to illustrate the effect of structural changes.High-level computations of a reaction path are restricted to a few photochemically relevant excited states that are involved in the main electronic transitions.The here presented simulations can complement the picture by employing a state manifold for the initial vertical photoabsorption.
The network for such a computation is displayed in Figure 2. Here, the Read External plugin is employed to read-in the results of an ORCA calculation in order to carry out a state-selective excitation and to visualize the resulting manyelectron dynamics.The dynamics is treated by the hybrid time-dependent density functional theory/configuration interaction (TDDFT/CI) formalism 40,68 : In a first step, the structure of guanine is optimized at CAM-B3LYP/cc-pVTZ level of theory using the ORCA package.From this geometry, a total of 10, 109 singlet excited states are obtained using LR-TDDFT in the Tamm-Dancoff approximation.The output, namely the molecular geometry, basis set definition, orbitals, pseudo-CIS eigenvectors, and the corresponding state energies is then read in and stored by the Read External module.
Once the wave function is available, the calculation of dipole matrices, μ i,j;q ¼ Ψ CIS i jb μ q jΨ CIS j D E , is done using the Dipole Matrix module.It provides useful information on the state properties.The transition dipoles on the off-diagonal lay out which states can be addressed by optical excitation.The permanent dipoles on its diagonal allow a first estimation of the character of an induced transition.The picture is often supplemented by a visual inspection of the involved orbitals, in Jellyfish done by the Plot Orbtials and States module.Details on how fast a state-specific transition can be obtained within a large excited states manifold and what pulse intensities or shapes are necessary, can only be answered by explicit time-dependent numerical simulations.In Jellyfish, they are available through the propagation module.In the following, an exemplary excitation to a specific state is outlined.
The first excited state of our LR-TDDFT calculation was selected as the target for a state-to-state transition.It is a bright state with reasonably large transition dipole moments: μ 0,1 ¼ À0:788, À 0:936,0:010 ð Þ ea 0 , thus making it easy to be addressed by linearly polarized laser field.
The target state is located 0:200 E h (5.442 eV) above the ground state and consists primarily of a HOMO-LUMO transition (82%).These two frontier orbitals are depicted Figure 4a,b.Excitation to this state involves a reorganization of the electronic π system with respect to the ground state as depicted by the difference density, Figure 4c.It reveals the overall rearrangement of π orbitals on almost all atomic centers upon excitation.The red isosurfaces show the withdrawal of electronic density while the green ones show the gain and, thus, the excitation can clearly be characterized as 1 ππ Ã ð Þ state.
To obtain this rearrangement of stationary densities, a full population inversion between ground and target state has to be induced by a laser pulse.This may be achieved by more or less complex shaped pulses in the following, however, a resonant excitation by a so-called π pulse with optimized parameters to obtain the targeted transition, is presented.
In the following, the photoexcitation by a cosine squared laser pulse with the duration of 2σ ¼ 2000ℏ=E h (48 fs) is presented.Regarding the conditions for a π pulse, f 0 ¼ πℏ=σ j μ i,f j, this duration results in a moderate maximum field amplitude of 0:002E h = ea 0 .The propagation was performed in 5000 time steps with a step size of 0:4ℏ=E h in the basis of the 1500 energetically lowest pseudo-CIS states.The time-dependent state populations, ij, and the time-dependent dipole moment as a result of the π pulse excitation are shown in Figure 5.Note in passing that for this figure, both the state populations as well as the time-dependent dipole moment have been exported from Jellyfish and plotted by an external program.As seen, the populations smoothly interchange with minute oscillations originating from the half-cycles of the oscillation laser field.At the end of the pulse, however, the target state population (orange curve) nearly reaches the maximum value of 1 and the ground state population (blue curve) nearly completely drops to 0. That means, the state-to-state transition is incomplete, which can be due to several effects.For one, the coupling of the permanent dipole moments to an intense external field can modify the transition energy so that it becomes time-dependent (dynamic Stark effect).Furthermore, a short laser pulse has a relatively broad spectral width so that also neighboring states can be reached by the pulse.The incompleteness in this particular case, however, is mostly due to multi-photon transitions to higher excited states, as seen in the sum population of ground and target state (green curve).Inspection of the entire state population revealed that a multitude of states at integer multiples of the incident laser frequency have gained very small occupations < 1% ð Þduring the pulse excitation.The time-dependent dipole moment, plotted in gray, starts at the total permanent dipole moment of the ground state (j μ 0 j¼ 2:549 ea 0 ) and shows distinct oscillations during the excitation.These oscillations coincide with the frequency and the evolving intensity of the laser pulse.At the end of the pulse (2000 a.u.) the time-dependent dipole moment has almost reached the value of the permanent dipole moment of the target excited state of j μ 1 j¼ 2:178 ea 0 .Since the excitation is not complete, by the end of the propagation an electronic wave packet has been created so that the dipole moment is noticeably oscillating (and would continue to oscillate even after the pulse is off ).
Aside from these quantifiable observables, the depiction of densities gives detailed information on the electronic redistribution in terms of their spacial extent and their localization at atomic centers.Especially the view of the evolution in time can be insightful.For the purposes of the program presentation, a short Jellyfish animation is provided.See Ref. 7 for a more comprehensive description and Jellyfish analysis of an excitation to a charge transfer state.
A video of four different density representations over the period of the central laser cycle from 900ℏ=E h to 1100ℏ=E h can be found online.This animation is also summarized in a sequence of snapshots in Figure 6 at significant points of the laser cycle.For a carrier frequency of ℏω ¼ 0:2 E h the zero crossings of the center optical cycle are roughly at t ¼ 984ℏ=E h , t ¼ 1000ℏ=E h , and t ¼ 1016 ℏ=E h , so that the minimum and maximum field strengths are reached around t ¼ 992ℏ=E h and t ¼ 1008ℏ=E h , respectively.
In the first column of the figure, the time-dependent difference density is shown where the red isosurface corresponds to a loss and the green isosurface to a gain of local density.The second column displays the gradient density with blue isosurfaces indicating a gain and pink a loss of electron density.In the last two columns, the NTO densities for hole (gray) and particle (blue) are rendered.
Before looking at the progressions on certain atomic centers, the four displayed densities seem to be divided into two groups.Apparent at first sight is that both NTO densities on the right, ρ NTO e=h , evolve much less dynamically but rather show a smooth and steady increase in volume, while the two kinds of difference densities on the left, Δ 0 and Δ t , F I G U R E 5 Time-dependent dipole moment (gray) and the population of the ground (blue) and target (orange) states of guanine during a state-to-state transition by a π pulse.
show more complex structures along rapid oscillations.Note, the stationary electron-hole NTO density and the difference density are in the CIS formalism mathematically equivalent.This is not the case during a laser excitation.The time-dependent NTOs show the formation of hole and particle in terms of the rising population in the excited states, cf., Figure 6 Thus, the NTOs are essentially excited state quantities/properties, since they are calculated with the ground state contributions are being removed to ensure orthogonality during time propagations (cf., Ref. 55).
However, the two types of difference densities are calculated directly from the field-dependent densities and therefore show pronounced oscillations that coincide with the slopes of the external field.The oscillations of Δ 0 mirror the oscillations of the dipole moments, compare Figure 6, while the induced flow of current at this time is visualized by the gradient density Δ t .
For the difference densities, the picture of p orbitals on atomic centers being perpendicular to the ring plane is true for certain times only.That is, when the slope of the optical cycles is zero, that is, when a laser cycle reaches its maximum or minimum.In between these significant points in time, the isosurfaces reflecting loss and gain have complex shapes similar to heavily distorted d orbitals with contributions of loss and gain at opposite sides of a single atomic center.
The participation of certain atomic centers in the density rearrangement is not equally distributed.For example, among the hetero atoms, the oxygen atom of the carbonyl group is showing more loss than gain on average, while the adjacent carbon atom gets involved only after the first third of the propagation (not shown) showing mostly the deposition of electronic density once the laser strength has become more intense.The same applies for the nitrogen atom N 7 in the 5-membered ring, it seems rather inert during the first third of the pulse.Furthermore, while the three hydrogen atoms of the five-and six-membered rings do not participate much, the two hydrogen atoms of the amino group show an active involvement in the rearrangement of electronic density.This activity during the pulse excitation is not apparent from the analysis of stationary densities (only slightly suggested by the HOMO and LUMO figures).

| Laser ionization in quantum dots
Dynamics simulations by means of propagating an electronic wave function can be done in explicitly atomistic descriptions but also by model potentials.To illustrate the capabilities of the Jellyfish package we present a model potential description of a QD and its light-induced ionization process.QDs are semiconductor nanomaterials.Their optical and electronic properties can be controlled through their size and shape or by gating through external electrodes and can be captured by effective mass model potentials, for example, Gaussian functions 54,69,70 of the form where D is the depth and b is a measure of the width of the QD in atomic units.Here the parameters for b and D are 0:2a À2 0 , 2 E h , respectively, which is supposed to resemble a gallium-arsenide QD from experimental data. 69Standard basis sets have been optimized for bound states of atoms of the periodic table but not for Gauss potentials.However, constructing custom-made basis sets following the systemics known for atoms can be done in the form of an eventempered basis set generated with the corresponding module in Jellyfish.Even-tempered basis functions centered at the coordinate origin are defined according to where a j is the angular momentum in the corresponding spatial direction.They cover the Hilbert space evenly, which allows them to represent both bound and continuum states.Due to their squared exponential form, they can be computed efficiently with the same integrators as for standard basis sets.For our example, two sets of even-tempered basis functions with total angular momentum a ¼ a x þ a y þ a z = 0, α = 0:07, β 0 = 1:5, and N = 11, as well as with angular momentum a = 1, α = 10 À4 , β 0 = 1:8, and N = 21 were chosen.Smaller values for β 0 , that is, more diffuse Gaussians, can lead to linear dependencies.Here, the model potential is considered for two electrons occupying both the single bound state available.From a full-CI calculation, a ground state energy of À1:149 E h was obtained and an ionization energy of ε HOMO ¼ 0:277 E h was determined using Koopmans' theorem.To simulate ionization in this system, a cos 2 laser pulse (Equation ( 1)) as shown in the bottom panel of Figure 7c was used with a frequency of ℏω ¼ 0:057 E h , corresponding to a typical 800 nm Ti-sapphire laser, with a maximum amplitude of f 0 ¼ 0:2 E h =ea 0 (1:4 Â 10 15 W/cm2) and a pulse duration of 700 ℏ=E h (17 fs).The laser intensity was chosen to promote an electron above the ionization barrier.Free propagation of an unbound electron, however, is restricted due to the use of atom-centered basis sets.To compensate for such numerical issues a spherical CAP in real space, cf., Section 3.1 was introduced to the total Hamiltonian.The CAP has a quadratic rise sitting on a distance from the origin of 15 a 0 .It has a quadratic fall-off until the maximum value of 10 E h is reached. 71The necessary integrals between the CAP and the basis functions were obtained by numerical integration on an even-spaced rectangular grid between À25a 0 ;25a 0 ½ with 200 grid points.Propagation was performed for 1600 time steps with a step size of Δt ¼ 0:5 ℏ=E h using all 5041 CI states as the basis.All steps for the calculations the evaluation were also performed directly in the JellyfishGUI.Here, the time-dependent populations were calculated and displayed as shown in the screenshot in Figure 7a.
States that have an overlap of at least 1% at any point in time to the total time-dependent wave function are displayed in the corresponding window.The highest of these states (state 76) has an excitation energy of 0:770 E h .All states except the ground state are continuum states and the ionization progression is reflected in the time-dependent F I G U R E 6 Sequence of density snapshots shows the excitation process of guanine by a resonant π pulse during the most intense central laser cycle.The first column shows the difference density with respect to the ground state, the second the difference with respect to the previous time step's density, the third the NTO hole density, and the fourth the NTO particle density.
norm, as displayed in the screenshot in Figure 3b.From the populations, it can be seen that already at a slow increase of the laser amplitude starting at about 100 ℏ=E h first transitions into energetically low continuum states (58, light green and 67, yellow) takes place.However, the ground state and these continuum states have very little overlap with the CAP, so that the norm barely changes.With higher amplitudes, starting at 200 ℏ=E h the electron is excited into energetically higher continuum states, which have a larger overlap with the CAP, that is, the continuum electron has a higher kinetic energy and the electron density of the corresponding states is thus removed from the system more quickly.With decreasing laser amplitude excitation into continuum states stops.At the same time, the electron density of previously populated continuum states is completely removed so that after the termination of the pulse only population in the ground state remains.

| High-harmonic generation spectra of H 2
Attosecond laser pulses, as used in the previous examples, utilize the nonlinear optical properties of materials and are generated by the process of high-harmonic generation (HHG).Since this also is primarily an electron dynamics process, the TD-CI method can used for its calculation.In the following the computation of molecular HHG spectra with Jellyfish is demonstrated on the example on H 2 .
Since HHG spectra of H 2 and the influence of different simulation parameters have already been studied theoretically in great detail, 2,72,73 to list only a few works, the following will focus on two possibilities for the calculation of this type of spectra utilizing the Jellyfish program.
Adopting the parameters used by White et al. 72 the H 2 bond length is set to R HH ¼ 1:4 a 0 and the Dunning basis set d-aug-cc-pVTZ with additional diffuse functions was used.In all calculations a cos 2 shaped laser polarized along the H 2 bond with a frequency of ℏω ¼ 0:057 E h , corresponding to an 800 nm laser, with a maximum amplitude of f ¼ 0:0534 E h =ea 0 (10 14 W/cm 2 ) and 20 optical cycles was used.This leads to a pulse duration of 2204:6 ℏ=E h (53:33 fs), the total propagation time, however, was 3000 ℏ=E h (72:57 fs) with a step size of Δt ¼ 0:01 ℏ=E h for the subsequent numerical Fourier transformation.
The HHG power spectrum P ω ð Þ for a specific polarization direction can be calculated from the propagation using the Fourier transform of the time-dependent dipole moment μ !t ð Þ as where t i and t f correspond to the initial and the final time.W is an optional window function that can help to minimize the noise in the spectrum.For this purpose, the Jellyfish package offers the option to employ a Hann window function.Furthermore, instead of the dipole moment, also the second derivative of the dipole moment with respect to time (dipole acceleration) is often used, which can also be selected in Jellyfish.This FT with options can be done in the same analysis window as shown in Figure 7 in the corresponding window tab.In Figure 8, three HHG spectra for molecular hydrogen are shown.The spectrum in blue was obtained straightforwardly from simulations with the molecular Hamiltonian only, while the spectra in orange and green are the result of additionally employing two types of CAPs.As already presented in the QDs example a real-space CAP can be included in the electron dynamics simulations (Section 4.2).][75] The second CAP implemented in Jellyfish, a CAP in energy space similar to the heuristic CAP presented by Klinkusch et al. 76 can also be used for HHG spectra.As in the ab initio version of Coccia et al., 77 each molecular orbital ϕ r ! with positive orbital energy ε p > 0 is assigned an inverse lifetime γ p as ε p ! ε p À iγ p =2.In the second quantization, this CAP operator can thus be written as As in the heuristic CAP, d is an escape length which in this example was chosen as d ¼ 20 a 0 .The calculated HHG spectra without and with spatial and energy CAP are shown in Figure 8.
Molecular HHG spectra reveal the linear and nonlinear response properties of the molecules exposed to an external light source.Usually, there is a sharp decrease of signal intensities after the fundamental band (first harmonic, linear response) followed by a plateau region of harmonics with similar intensities until a cutoff is reached where the signal strengths steadily decrease.This progression is often explained by the semi-classical three-step model of Corkum. 78The harmonic signals themself are the result of intense optical laser cycles.Upon ionization, the electronic wave packet is accelerated by the rising slope of a laser half-cycle (step one) followed by an acceleration back to the parent ion by the second half-cycle of the opposite sign (step two).The gain of kinetic energy of the wave packet is released upon recollision with the parent ion (step three) in the form of multiple integers of the incident laser frequency.The outlined intensity progressions of the harmonic signals can be seen in all three spectra shown in Figure 8.After the fundamental, a series of odd harmonics is identified-as a consequence of the selection rules for optical excitations, only odd harmonics can be obtained in a homonuclear diatomic molecule.The spectrum in blue (where no CAP was used) appears much noisier as it also records a multitude of electronic transitions to and within highly excited states resulting in sharp bands at high energies in addition to the harmonic signals.The presence of discrete highly excited states above the ionization potential is an artifact of the usage of finite atom-centered Gaussian basis sets that have been optimized for the description of atomic and molecular bound states.Hence, the inclusion of an absorbing potential for such high-lying excited states is necessary when strong and intense laser fields are applied.
The effect of the energetic and the real-space CAP is well seen in the spectra plotted in orange and green.The harmonic signals are much more clearly identifiable since excitations to states in the continuum are suppressed, even differences in their shape widths can be analyzed.However, up to the 15th harmonic, both types of absorbers deliver almost identical signal intensities.For the following harmonics, the signals obtained by the spatial CAP are less strong.Furthermore, the green spectrum shows additional signals at positions that do not fit the odd-harmonics frequencies.Here, the difference between the two absorbing potentials to capture an ionization process is manifesting.For the HHG spectra, the energetic CAP has the advantage of assigning lifetime to all states above the ionization potential.In the case of the real-space CAP in the form of a sphere around the molecule, not all high-lying excited states may have an overlap with the absorber.This is again a well-known artifact of the finite size of the atom-centered basis set that creates not only energetically, but also spatially localized excited states within the ionization continuum.As a consequence, in the case of the spatial CAP, population can be "trapped" in such localized states and the interaction with the laser field causes additional false signals.Further adjustments to the basis set and inclusion of more diffuse functions may improve the HHG spectra quality but at the expense of calculation time.However, spatial CAPs have the advantage of better capturing angle-dependent ionization yields as shown by Schlegel and coworkers. 75,79

| CONCLUSION
Electron dynamics is a very active research field where experimental techniques and theoretical developments are rapidly evolving.Unlike in the related field of electronic structure methods, the existing dynamics programs are mostly limited to a few specific methods and systems.With Jellyfish, we have developed a program that provides a flexible, user-, and developer-friendly environment.This is possible on the one hand by the graphical user interface which F I G U R E 8 HHG spectra of H 2 at a laser intensity of 10 14 W/cm 2 calculated with TD-CIS in the d-aug-cc-pVTZ basis set without an ionization model (blue) and with two different CAP operators to account for ionization process.

F I G U R E 1
Program architecture of Jellyfish.JellyfishCore contains the basic definitions of the plugin and module system and is linked to both plugins and programs to enable the necessary communication.The plugins contain the individual modules for Jellyfish and are dynamically loaded into the programs as a shared library.F I G U R E 2 Screenshot of the main window of the JellyfishGUI.Shown is the network to read in the electronic structure calculation from an ORCA log file to carry out a laser-driven electron dynamics simulation: (a) NodeEditor, central window in which modules can be selected, connected, disconnected, and executed; (b) List of available plugins and their modules; (c) Detailed settings for the modules chosen; (d) Project log file; (e) Tabs to open and switch between projects.

F I G U R E 3
Screenshot of the visualization module, Plot TDCI Densities, for time-dependent electron densities.(a) In the central window the wave function objects are rendered.In this case, the difference density, Δ 0 , between the ground state and the time-dependent density is shown; (b) The left panel allows the selection of density and the time frame.This selection menu changes depending on the chosen type of density; (c) Control options for the visualization, like the color and number of isosurfaces, as well as graphics settings controlling the picture quality and export functions; (d) The lower frame displays the log file of executed actions and error messages.

F I G U R E 4
Selected static properties of the first excited state of the guanine molecule: HOMO, LUMO, and difference density between ground and first excited state, Δ 10 .

F I G U R E 7
JellyfishGUI screenshots from subwindows for the evaluation of driving laser field and time-dependent properties.(a) Population of all states with an overlap of at least 1% with the wave function, (b) norm, and (c) laser pulse.