Szczepanik Research Group

Department of Theoretical Chemistry
Faculty of Chemistry, Jagiellonian University
Gronostajowa 2, 30-387 Krakow, Poland
Tel: (+48) 12 686 23 90
E-mail: dariusz.szczepanik@uj.edu.pl

HOME      PROJECTS      PUBLICATIONS      runEDDB      molANO      X

About EDDB/BDF       Download       Tutorials       Manual


◼ Manual (2026-Jun-04)

Table of Contents
1. Introduction
2. Generating Input Files
3. Analysis Modes
4. Command-Line Reference
5. Output Files
6. Differential Maps
7. Understanding the Output
8. Known Issues and Limitations
9. Performance Tips
10. Citations
11. Contact

1. Introduction

runEDDB implements the Electron Density of Delocalized Bonds (EDDB) method and the corresponding electron-density partitioning schemes. The present release is a complete rewrite in Fortran 2008 of the original runEDDB R-script. Thanks to a modern linear-algebra backbone (BLAS-3 DGEMM, DSYGVD, DSYEVD, DPOSV, CGS2, etc.), the program now processes systems containing thousands of atoms on a single workstation in seconds to minutes. It is distributed as ready-to-run, pre-compiled binaries for all three major desktop platforms — Linux (x86_64), Windows (x86_64), and macOS (arm64). See the Download page for the current release and per-platform installation instructions.

The Bond-Order Projection (BOP) algorithm at the heart of EDDB requires the input wavefunction to be expressed in a basis of compact, valence-localized atomic orbitals with a hierarchical, L-resolved subshell orthogonalization — the Natural Atomic Orbital (NAO) basis. The NAO transformation maps the underlying Gaussian basis onto an atomic-orbital basis with a clean core (Cor) / valence (Val) / Rydberg (Ryd) hierarchy. Several implementations of the original definition are in routine use today (e.g., NBO7, JaNPA); runEDDB ships with its own internal implementation alongside two compact subsets — the Natural Minimal Basis (NMB = Cor ∪ Val) and the program-default Natural Valence Basis (NVB = Val).

2. Generating Input Files

The current version of runEDDB directly reads seven input formats — six wavefunction formats produced by external quantum-chemistry programs and one geometry-only XYZ format processed by the built-in GFN2-xTB engine:

FCHK (.fchk)

The most common input format. Supported routes:

Most of our FCHK testing focused on Gaussian and Q-Chem outputs. For post-HF methods the relaxed correlated density block has to be requested explicitly — in Gaussian via density=mp2 (MP2) or density=current (CCSD, CI, CASSCF) on the route line; in Q-Chem via cc_ref_prop = true for the CC family. (Q-Chem does not currently emit the MP2 density to FCHK; for MP2 use Gaussian or ORCA.)

Molden (.molden)

Molden is one of the most popular formats for molecular density and orbital visualization. Many quantum-chemistry programs can produce it through different mechanisms — please consult the manual of your particular package for the exact keyword. Multiwfn is also a convenient universal converter to and from Molden, and runEDDB has been validated against Multiwfn-generated files.

ORCA (.json or .molden)

For ORCA we strongly recommend using .json rather than .molden. The JSON format produced by ORCA carries a more complete description of the wavefunction — including ECP definitions and other metadata — that the Molden specification cannot represent. We have validated runEDDB against both formats, generated by ORCA's standard orca_2mkl and orca_2json post-processors.

Turbomole native (.tmol)

For Turbomole calculations, runEDDB reads a single .tmol bundle that you assemble by concatenating the four key files from the Turbomole working directory:

 terminal
$ cat control coord basis mos > molecule.tmol

The parser handles RHF, UHF, ROHF, and DFT wavefunctions in Turbomole's spherical basis representation, with angular momenta up to k (L = 7) — covering def2-{S,T,Q}ZVPP and cc-pVDZ through cc-pV8Z. The native Turbomole $ecp block is parsed automatically — ECP-bearing transition-metal calculations load directly without conversion. Cartesian d/f basis storage (forced by $pople CAO or $cao) is not supported — use Turbomole's default spherical mode.

C1 symmetry is required for the .tmol lane. Under any non-trivial point group Turbomole writes the wavefunction in two ways that need irrep-aware bookkeeping not present in the .tmol parser:

The fix is simple: re-run the SCF with $symmetry c1 in the control file (or remove the keyword). Population analysis with EDDB doesn't benefit from point-group symmetry being on during the SCF write-out, so this is the standard workflow. If a non-C1 .tmol is passed, the parser stops with an explanatory ERROR; the run does not silently produce wrong densities. The other input formats (.fchk, .json, .molden, .49, .gms) write full-AO MOs and global orbital ordering regardless of point group, so they have no equivalent restriction.

GAMESS-US (.gms)

GAMESS-US output files (the standard text log produced by rungms) are read directly. Just rename or symlink the log to a .gms extension — PUNCH (.dat) files are not used. Supported run types: RUNTYP=ENERGY, OPTIMIZE (final equilibrium geometry is auto-located), and SADPOINT (transition state). Supported wavefunction types: RHF, UHF, ROHF, and MCSCF natural orbitals. Spherical-harmonic basis (ISPHER=1) is recommended; standard Cartesian basis is also supported. The default GAMESS print precision (F11.6) is sufficient; for high-precision basis-set studies, set DGRID=.T. in $CONTRL for F17.10 MO output.

GAMESS $ECP blocks are parsed automatically (handles both WITH ZCORE n AND LMAX l and the shared-definition shorthand ARE THE SAME AS ATOM N). Symmetry-truncated basis sections (any GAMESS run with COORD=UNIQUE under Cs, C2v, D2h, etc.) are reconstructed automatically — the parser detects the symmetry-implied gaps in the printed basis section and replicates each per-atom block for its unprinted symmetry-equivalent partners.

NBO archive (.49) — legacy compatibility

NBO7 .49 archives are still accepted in NAO-only mode (population analyses + BOP only; cube grids and orbital FCHK exports are disabled because NBO uses a different AO-basis convention than runEDDB's native pipeline). The .49 path is provided for backward compatibility and will be removed in a future release once the native NAO driver completes its acceptance suite. Users are strongly encouraged to switch to .fchk / .json / .molden / .tmol / .gms as soon as practical.

XYZ geometry (.xyz) — built-in GFN2-xTB engine

One of runEDDB’s most distinctive capabilities is that it can analyse chemical bonding straight from a bare set of atomic coordinates — for example a structure solved by X-ray diffraction, a molecular-dynamics snapshot, or a hand-built model — with no separate quantum-chemistry calculation. Give it a plain .xyz geometry and it builds the electronic structure itself, on the fly, using a fast built-in GFN2-xTB engine (a semi-empirical tight-binding model; Bannwarth, Ehlert & Grimme, JCTC 15 (2019) 1652), and then runs the complete EDDB analysis on top of it.

The practical payoff: in a single command, and in seconds rather than the hours a full DFT study would take, you obtain a map of where the electrons delocalize — a direct picture of the molecule’s resonance / conjugation pattern and bonding skeleton — without any external software, wavefunction file, or DFT step. This makes runEDDB a quick, stand-alone way to look at the bonding in a freshly determined or modelled structure; promising findings can then be confirmed with a higher-level wavefunction (DFT, post-HF, …) run through one of the other input lanes — again analysed by runEDDB — for a fully validated picture. The entire XYZ pipeline (GFN2-xTB SCC in its minimal valence vSTO-nG basis → AO→NAO localisation → BOP and EDDB analysis) runs end-to-end on the plain XYZ file, and the internal engine is validated against the reference tblite implementation across 81 systems covering 26 elements (H…W) at chemical accuracy.

About the GFN2-xTB basis (vSTO-nG). GFN2-xTB is a minimal valence tight-binding method: its atomic-orbital basis has just one Slater-type orbital per valence shell (e.g. one s and one p for C, N, O; one s for H; with a d shell added for heavier elements), and each Slater orbital is expanded as a short contraction of n Gaussians (an STO-nG fit). We abbreviate this vSTO-nG = valence STO-nG. It is therefore a compact, single-zeta-quality valence basis — not a double-zeta or polarized basis — which is exactly why an NAO/NBO-style decomposition of a GFN2 wavefunction looks minimal-basis-like. This keeps GFN2-xTB fast and is appropriate for the qualitative bonding screening the .xyz lane is meant for; for quantitative work, validate with a larger-basis DFT / wavefunction calculation through one of the other input lanes.

The XYZ file uses the standard XYZ format, with a runEDDB-specific optional extension on line 1 that carries the molecular charge and spin multiplicity:

 molecule.xyz
N_atoms [ charge multiplicity ]
comment line
Sym1   x1   y1   z1
Sym2   x2   y2   z2
...

The optional charge and multiplicity tokens on line 1 (after the atom count) default to 0 and 1 respectively, so plain XYZ files (with only N_atoms on the first line) are processed as neutral closed-shell singlets by default. Coordinates are in Ångström. Total energies reported in the verbose output of the GFN2-xTB stage match Grimme’s reference xtb executable to within ~10−4 Eh, allowing direct sanity-checking against an external xtb run on the same geometry.

Trajectories and geometry alignment (--xyz-*)

A .xyz file may hold a trajectory — many structures concatenated one after another (e.g. an MD run or an optimisation path). By default runEDDB analyses the first structure and prints a note reporting how many it found; use --xyz-frame N to pick the N-th structure instead (1-based). Requesting a frame that does not exist is a clear error reporting the actual count.

Two further options rigid-translate the input geometry (bond lengths, angles and dihedrals are untouched, so the EDDB analysis is unchanged — only the absolute coordinates written to any exported cube / FCHK move). They apply to every input format, not just .xyz (legacy .49 excepted):

Trajectory tip — comparable cubes across frames. Centering each frame independently with --xyz-center would re-centre every structure on its own (drifting) centroid, so the resulting cubes would not line up. Instead, centre one reference frame with --xyz-center geom -v, read off the printed translation vector, then run every frame with that same fixed vector via --xyz-translate (optionally combined with --cube-size 0 reference.cube to clone the grid) — the frames then share one coordinate frame and the EDDB cubes are directly comparable / subtractable.

Both closed- and open-shell systems are supported. Closed-shell molecules (the default, multiplicity = 1) are treated as restricted; any higher multiplicity is solved as a true spin-polarised unrestricted (UHF) calculation (spin-polarized GFN2-xTB / spGFN2-xTB; see the method references in §10) — just put the spin multiplicity on line 1 (e.g. 2 for a doublet radical, 3 for a triplet, 4 for a quartet, …).

For any open-shell run, runEDDB prints a spin-density topology table right after the GFN2-xTB step: the per-atom spin population, i.e. where the unpaired electrons actually sit. This is a quick way to read off the magnetic / radical centres of a structure — and a practical screen for diradical / broken-symmetry character: run the molecule as a triplet and, if the spin localises on two well-separated atoms or fragments, that is a strong hint of an open-shell singlet (diradical) ground state worth investigating. The --output-no export adds complementary diradical / multireference diagnostics — including the ⟨S²⟩ spin-contamination value, which counts roughly how many electron pairs are genuinely unpaired.

Diradicals / broken-symmetry singlets. GFN2-xTB’s simplified spin treatment cannot reliably converge a genuine broken-symmetry singlet — it tends to collapse back to the closed-shell solution or stop at partial open-shell character. Use the triplet spin-topology read-out above as a qualitative pointer only, and confirm any suspected diradical / broken-symmetry singlet with a proper broken-symmetry DFT calculation (e.g. CAM-B3LYP), which can then be analysed by runEDDB through one of the wavefunction input lanes.

Note — the GFN2-xTB lane is a qualitative screening tool. GFN2-xTB does not capture long-range exchange correctly, so the resulting EDDB indices (NOBD occupancies, CDI, EDDBG/F/P populations) are systematically inflated by ~10–20% relative to range-separated hybrid DFT. The lane is excellent for surveying the qualitative bonding topology — identifying which atoms / orbitals participate in delocalization, which interactions look chemically interesting, and which pathways are worth a closer look — on systems for which no wavefunction file is yet available (X-ray geometries, MD snapshots, hand-built coordinates). Chemically interesting findings should always be confirmed with a follow-up ab initio (DFT, post-HF, …) single-point on the same geometry, processed through one of the other input lanes. See the xTB tutorial for a worked example comparing GFN2-xTB and CAM-B3LYP/def2-SVP on the same metallacycle.

f-block elements. GFN2-xTB’s spin-polarisation parameters are defined only for s, p and d shells, so the open-shell character of lanthanide / actinide f electrons is not described — this is a property of the GFN2 model itself, shared with the reference tblite. f-block compounds still run, but any spin / delocalization assigned to the f shell should be interpreted with care, and a few lanthanide halides cannot even be cross-checked against tblite because the reference implementation does not handle them.

MDL Molfile (.mol) and SDF (.sdf)

The MDL Molfile (.mol) and its multi-record container SDF (.sdf) are a second geometry-only lane: like .xyz they carry just a structure (and, optionally, formal charges / radical flags), so runEDDB builds the electronic structure with the same built-in GFN2-xTB engine and runs the full EDDB analysis on top — everything said above about the GFN2-xTB lane (the vSTO-nG basis, the qualitative-screening caveat, open-shell handling, the geometry options below) applies unchanged. .mol / .sdf are ubiquitous in cheminformatics (PubChem, ChEMBL, RDKit, drawing programs), so this lets you analyse a structure straight from those sources.

runEDDB reads the standard V2000 connection table. Coordinates are in Ångström; the bond block is ignored (GFN2-xTB rebuilds connectivity itself). Molecular charge and spin are taken from the property block: total charge is summed from the M  CHG entries (falling back to the legacy per-atom charge column if absent), and multiplicity from M  RAD radical flags (doublet → 1, triplet → 2 unpaired electrons); when no radical flag is present the multiplicity follows the electron-count parity (even → singlet, odd → doublet). Two cases are refused with a clear message: the free-form V3000 format (re-export as V2000), and files explicitly flagged 2D (the dimensional code on line 2) — GFN2-xTB needs real 3D coordinates.

An .sdf file may hold many structures (conformer sets, compound libraries), each record separated by a $$$$ line. As with multi-frame .xyz (see above), runEDDB analyses the first record by default and prints a note with the count; --xyz-frame N selects the N-th record. The --xyz-center / --xyz-translate geometry options apply here too.

Tripos MOL2 (.mol2)

The Tripos .mol2 format is another geometry-only lane feeding the same GFN2-xTB engine. runEDDB reads the @<TRIPOS>MOLECULE atom count and the @<TRIPOS>ATOM records; each element is taken from the SYBYL atom type (C.3 → C, N.ar → N, or a bare symbol such as Li), and the bond block is ignored. Coordinates are in Ångström. Multiple @<TRIPOS>MOLECULE records are treated as a multi-structure file (--xyz-frame N selects one). MOL2 carries only per-atom partial charges, not a formal total charge or spin, so runEDDB defaults to a neutral closed-shell singlet; set a different state with --xyz-charge / --xyz-mult (below).

Protein Data Bank (.pdb / .ent)

The PDB lane reads ATOM / HETATM records (coordinates in Ångström). The element is taken from the dedicated element columns (77–78) when present — as written by Open Babel, PyMOL, etc. — and otherwise guessed from the atom name. Only the first alternate location (altLoc blank or A) of any disordered atom is kept. A multi-MODEL file (NMR ensembles, MD frames) is multi-structure: the first model is used by default, --xyz-frame N selects another. PDB carries no charge or spin, so the default is a neutral closed-shell singlet (override with --xyz-charge / --xyz-mult).

Hydrogens. Crystallographic PDB structures often omit hydrogens. GFN2-xTB needs a complete all-atom geometry, so a carbon-containing PDB with no hydrogens triggers a warning — add hydrogens first (e.g. obabel in.pdb -O out.pdb -h, or in PyMOL / Avogadro). Pure inorganic / metal structures that legitimately have no hydrogens are not flagged.

Overriding charge and multiplicity (--xyz-charge / --xyz-mult)

For any geometry input (.xyz, .mol, .sdf, .mol2, .pdb) the molecular charge and spin multiplicity can be set on the command line with --xyz-charge N (a signed integer) and --xyz-mult M (2S+1, a positive integer). When given they take precedence over whatever the file specifies (the .xyz header, or the M  CHG / M  RAD records of a molfile). For .mol2 and .pdb, which carry no usable charge / spin information (default: neutral singlet), these flags are the way to study an ion or an open-shell state — e.g. --xyz-charge -1 --xyz-mult 2 for a radical anion. They are rejected for the wavefunction lanes (.fchk / .json / .molden / .tmol / .gms / .49), where the charge and spin come from the wavefunction itself.

Crystallographic data (.cif) — convert to .xyz first

Crystallographic Information Files (.cif) are periodic data — a unit cell, an asymmetric unit (often only a fragment of a molecule, sitting on a symmetry element), and symmetry operators — not a ready molecular geometry. Turning that into a single, complete, chemically-sensible molecule (apply symmetry, grow the molecule across cell boundaries, resolve disorder and partial occupancies) is exactly what dedicated crystallography programs already do well. So rather than re-implement that, the recommended route is to export a complete molecule to .xyz from the tool you already use, then analyse that .xyz with runEDDB:

Once you have the molecule as .xyz, the full GFN2-xTB → EDDB pipeline runs on it directly (add hydrogens first if the structure lacks them, and set the charge / multiplicity with --xyz-charge / --xyz-mult if needed). This keeps the resonance / delocalization analysis one short step away from raw crystallographic data without duplicating the symmetry-and-disorder machinery your crystallography software already provides.

Effective Core Potentials (ECP)

Heavy elements are routinely treated with Effective Core Potentials (ECPs), which replace the chemically inert inner shells with an analytical pseudopotential. runEDDB consumes already-converged MOs and never evaluates ECP integrals itself; what each parser must do is read the per-atom effective nuclear charge Zeff = Zatomic − Ncore,replaced so that atomic charges, total molecular charge, and the per-atom population tables balance. The NAO classification (Core / Val / Rydberg) also uses each atom's Z to decide how many core/valence shells to expect, so getting Zeff right is essential — a wrong value silently misclassifies orbitals and can drop the heavy atom from the per-atom tables.

How each parser obtains Zeff:

Once Zeff is known, the NAO machinery automatically labels the surviving valence shells with the chemically correct principal quantum numbers. For example, tungsten under def2-ECP60 shows its valence as 6s, 5d1..5d5 (the actually-present W valence shells) rather than the misleading 1s, 3d1..3d5 labels that an ECP-blind labelling would produce.

Supported wavefunctions:

3. Analysis Modes

runEDDB offers three scope modes, two NOBD-dissection options, and three analysis-basis subsets.

Scope modes (mutually exclusive):

NOBD dissection (can combine with any scope):

When dissection is active, the atom table and summary show only the density of the selected NOBDs; the star marker (*) in the NOBD table flags them.

Analysis basis subset, selected by --basis VALUE (or -b VALUE). VALUE is one of:

VALUE is fully case-insensitive (Nao, NMB, nVb, etc. all work).

.49 input note. The NBO7 .49 archive flags every NAO as either NMB or Rydberg, with no further Cor/Val partition. -b nvb therefore cannot be honoured on .49 input and silently falls back to -b nmb. When -b nvb is requested explicitly, a one-line clarifier — (.49 archives carry no Cor / Val split) — is printed directly under the Target NAO representation: entry in the Stage 2 detail block, aligned with the basis-label value. The clarifier is verbose-only (-v), keeping the default-mode output clean. To get a real Cor/Val split, use the matching native input (FCHK / JSON / Molden).

4. Command-Line Reference

Basic usage:

 terminal
$ runEDDB --input FILE  [ Options ]

The input file (.fchk, .json, .molden, .tmol, .gms, .49, .xyz, .mol, .sdf, .mol2, or .pdb) is the only required argument; all other options have sensible defaults. The --help screen prints a compact, grouped summary (with the export options collapsed to --output-* and the engine controls to --xtb-*); the table below documents every option in full:

 command-line options
    --input, -i FILE          Input file (.fchk / .json / .molden / .tmol / .gms / .49 / .xyz / .mol / .sdf / .mol2 / .pdb)

    --xyz-frame N             Select structure N from a multi-frame .xyz/.sdf/.mol2/.pdb (default: 1)
    --xyz-center [geom|mass]  Recenter geometry at origin: geometric centroid (geom, default) or mass
    --xyz-translate x,y,z     Apply a fixed translation vector (Angstrom) to every atom
    --xyz-charge N            Override molecular charge      (geometry input only)
    --xyz-mult M              Override spin multiplicity 2S+1 (geometry input only)

    --global, -g              EDDB_G - entire molecular system (default)
    --fragment, -f ...        EDDB_F - molecular fragment, e.g., -f 1:6,10,12
    --pathway, -p ...         EDDB_P - specific bond pathway, e.g., -p 1-2-3-4-5-6-1

    --select-nobd, -s ...     Select NOBDs for dissection, e.g., -s 1:3,4,5
    --auto-pi, -a             Auto-select pi-type NOBDs to get pi-density component

    --basis, -b ...           NAO basis subset: nao | nmb | nvb (default: nvb; .49 input: nmb)

    --output-eddb ...         Export EDDB to Gaussian CUBE file
    --output-edlb ...         Export EDLB to Gaussian CUBE file
    --output-bdf  ...         Export BDF  to Gaussian CUBE file
    --output-ed   ...         Export ED   to Gaussian CUBE file

    Differential maps (df = f(this) - f(REF), point-by-point on REF's grid; not .49):
    --output-ded   REF [OUT]  Differential ED:   dED   = ED(this)   - ED(REF)   to CUBE
    --output-deddb REF [OUT]  Differential EDDB: dEDDB = EDDB(this) - EDDB(REF) to CUBE
    --output-dedlb REF [OUT]  Differential EDLB: dEDLB = EDLB(this) - EDLB(REF) to CUBE
    --output-dbdf  REF [OUT]  Differential BDF:  dBDF  = BDF(this)  - BDF(REF)  to CUBE

    --output-nobd ...         Export NOBDs to Gaussian FCHK file
    --output-no   ...         Export NOs of the spinless 1-density to Gaussian FCHK file
    --output-ao   ...         Export NAOs  to Gaussian FCHK file
    --output-xtb  ...         Export the GFN2-xTB wavefunction to Gaussian FCHK file (.xyz input)
    --output-moa  ...         Export per-MO ED/EDLB/EDDB populations to a text file (HF/DFT; not .49)

    GFN2-xTB engine (.xyz input only):
    --xtb-etemp K             Electronic (Fermi smearing) temperature in K (default ~300; 0 = off)
    --xtb-conv X              SCC convergence threshold on max|dq| (default: 1e-4)
    --xtb-maxiter N           Maximum number of SCC iterations (default: 1500)
    --xtb-noaes               Disable anisotropic electrostatics (faster; not EDDB-grade)
    --xtb-nospin              Open-shell: use ROHF instead of the default spin-polarised UHF

    --cube-size, -c ...       Grid points per axis, ~|N|^3 total (default: 100)
    --ncores,    -n ...       Number of CPU cores (default: 1)
    --quiet,     -q           True-quiet: banner + summary block + closing time-stamp only
    --verbose,   -v           Show all tables and full stage details
    --help,      -h           Print this help message

runEDDB has three verbosity tiers:

Independently of the verbosity tier, the per-NOBD table is also printed whenever -s, -a, or --output-nobd is in use, the per-orbital table whenever --output-ao is in use, and the NO-resolution table (plus the multireference diagnostics block under it) whenever --output-no is in use — in default mode. (In true-quiet mode no tables print regardless of those flags.)

GFN2-xTB engine controls (--xtb-*)

These flags tune the built-in GFN2-xTB SCC and apply only to .xyz input (they are ignored for the wavefunction-file lanes). The defaults are chosen to be robust and EDDB-appropriate, so most runs need none of them; they exist for the occasional stiff or unusual system.

If the SCC does not converge, runEDDB stops and prints a short reminder suggesting exactly these knobs (a higher --xtb-etemp, a larger --xtb-maxiter, or a looser --xtb-conv), so you rarely need to remember them in advance.

5. Output Files

Outputs are generated only when explicitly requested with --output-*.

Gaussian cube files (.cube) — 3D electron density on a PCA-aligned grid; viewable in VMD, Chemcraft, GaussView, Avogadro2, IQmol, or any cube-compatible program.

Grid quality is controlled by --cube-size N (default: 100), with the grid covering the molecular bounding box (or a fragment subset, when -f / -p is in effect) plus a small margin. The sign and zero of N select between three modes:

Gaussian FCHK files (.fchk) — viewable in GaussView, Avogadro2, IQmol, or Chemcraft.

For all three FCHK exports the per-orbital occupation (NAO population, NOBD occupation, or NO occupation number) is stored in the Alpha Orbital Energies field, pre-scaled so that the value Avogadro2 shows in its eV-labelled energy column reads back numerically as the unscaled occupation. So a doubly-occupied orbital appears as "2.00 eV", a singly-occupied one as "1.00 eV", and so on — the unit label is technically eV but the number is the actual occupation. Other viewers (GaussView, IQmol, Chemcraft) that read FCHK energies in hartree will show the raw small values; if you need the actual occupations there, multiply by 27.2114.

Text file (.txt) — MO-resolved analysis.

Filename rules. When no explicit filename is given, it is derived from the input. The orbital FCHK suffix matches the active basis subset (Section 3):

 filename examples
Input:   molecule.fchk        (or .json / .molden / .49)
ED:      molecule.ED.cube
EDDB:    molecule.EDDB.cube
EDLB:    molecule.EDLB.cube
BDF:     molecule.BDF.cube
NOBD:    molecule.NOBD.fchk
NO:      molecule.NO.fchk
NAO:     molecule.NVB.fchk    (-b nvb, default)
NAO:     molecule.NMB.fchk    (-b nmb)
NAO:     molecule.NAO.fchk    (-b nao)
xTB:     molecule_xtb.fchk    (--output-xtb, from molecule.xyz)
MOA:     molecule.MOA.txt     (--output-moa)
dED:     molecule.dED.cube     (--output-ded REF)
dEDDB:   molecule.dEDDB.cube   (--output-deddb REF)
dEDLB:   molecule.dEDLB.cube   (--output-dedlb REF)
dBDF:    molecule.dBDF.cube    (--output-dbdf REF)

The same default filenames are used regardless of --global / --fragment / --pathway scope and regardless of --auto-pi / --select-nobd dissection — the scope and dissection are reflected in the file contents and in the summary line labels (EDDB_G / EDDB_F / EDDB_P), not in the filename. Provide an explicit filename after any --output-* flag if you need a custom name.

6. Differential Maps (ΔED / ΔEDDB / ΔEDLB / ΔBDF)

A differential map subtracts one runEDDB density field from another, point-by-point on a shared grid, to reveal how the electron structure reorganizes between a reference and a perturbed system — a substituent change, a hydrogen / halogen bond, a π···cation or π···HF contact, solvation, or a redox step. The differential of the bond-delocalization field, ΔBDFπ, maps π-electron (de)localization changes directly in the language of resonance; the same operation applies to every density in the EDDB decomposition, so there is one differential flag per density:

In every case REF is a cube of the same density written by an earlier run: it supplies both the grid (origin / points / steps, cloned exactly) and the reference values, so the subtraction is point-for-point, no external cube-arithmetic tool (Multiwfn, VMD, …) is needed, and no separate --cube-size 0 is required. OUT is an optional output filename (defaults <input>.dED.cube, .dEDDB.cube, .dEDLB.cube, .dBDF.cube). The differential maps run on the native lane (not the legacy .49).

π-resolved (ΔBDF) versus full (ΔEDDB / ΔEDLB). ΔBDF is normally a π-only quantity, so it should be built from the π-delocalization (π-NOBDs; see the workflow below). When a clean π-dissection is awkward, ΔEDDB and ΔEDLB on the full delocalized / localized densities need no π-selection at all and still show where delocalization is gained or lost — frequently the easier and more robust route.

Recommended ΔBDFπ workflow. The map is only as good as the set of NOBDs included, and the BDF cubes for both systems must be built the same way. (1) compute and save the reference π-BDF cube, e.g. runEDDB -i ref.fchk --auto-pi --output-bdf ref.BDF.cube; (2) for the perturbed system run runEDDB -i pert.fchk --auto-pi --output-dbdf ref.BDF.cube pert.dBDF.cube. Capture the π-system consistently in both molecules: --auto-pi selects π-symmetry NOBDs automatically, but a perturbation can change how many NOBDs carry the π-delocalization (and the classifier may miss one in a given system). Inspect the NOBD-resolution table for each molecule (printed with --select-nobd / --auto-pi / --verbose) and, if needed, replace --auto-pi with an explicit --select-nobd LIST so that all NOBDs contributing to the π-(de)localization are included in both runs — otherwise ΔBDF shows a spurious feature wherever a real contributor was dropped from one side. Typical ΔBDFπ isosurfaces are around ±0.010.

Same coordinate frame, same scope. A point-by-point subtraction is only meaningful if the two systems share one frame, so align them before the wavefunction calculation — keep the common fragment’s Cartesian coordinates fixed (runEDDB translates, but does not rotate, a supplied wavefunction). runEDDB warns if the reference and current geometries do not overlap. Both runs must also use the same scope (--global / the same --fragment / --pathway) and the same density flavour, so that their grids, NOBD sets and densities are comparable.

Large, flexibly-reorganizing systems — focus on a fragment. When a perturbation reshapes a remote reactive site so the whole molecule no longer overlaps, restrict the analysis to the conserved region of interest with --fragment (or --pathway). This focuses the EDDB density, the cube grid and the geometry-overlap check on the selected atoms only, so the differential reports how that fragment’s π-system responds to a change happening elsewhere — without being penalised for motion in the part you are not studying. The overlap check is purely spatial and element-based, so the reference and perturbed structures may number their atoms differently: as long as the studied fragment is physically superimposed, the overlap is recognised. (Keep that fragment’s coordinates frozen between the two systems, so the differential reflects electronic reorganization rather than atomic displacement.)

7. Understanding the Output

runEDDB prints up to three population tables. By default only the per-atom table; -v prints all three; -q prints none (true-quiet suppresses tables, per-stage messages and timings — only the banner, the EDDB summary block, and the two bookend stage lines around it are emitted).

NOBD-resolution table: appears with --auto-pi, --select-nobd, --output-nobd, or --verbose. Each row is one NOBD (Natural Orbital for Bond Delocalization):

Orbital-resolution table: appears with --output-ao or --verbose. One row per NAO in the active subset; row index matches the MO index in the --output-ao FCHK 1:1, so the table doubles as a navigation key when inspecting orbitals in Avogadro2 / IQmol / GaussView.

NO-resolution table: appears with --output-no. One row per natural orbital with non-trivial occupation, sorted by descending NOON. Each row reports the NOON itself plus the breakdown of its occupation into Core (Cor), Valence (Val), and Rydberg (Ryd) NAO contributions, with percentages. Row index matches the MO index in the --output-no FCHK 1:1. Below the table, three single-line multireference-character indicators are printed: NUMayer = Σi ni(2−ni) (zero for a pure single Slater determinant; rises with correlation), NUHG = Head-Gordon's smoothed variant Σi ni2(2−ni)2 (less inflated by the small-correlation tail of large basis sets — the most-cited NOON-only multireference indicator), and the Yamaguchi y0 biradical index (0 for closed-shell single reference; 1 for an ideal singlet diradical; intermediate for partial diradical character). Useful as a quick MR-character sanity check on post-HF or broken-symmetry inputs.

Atomic-resolution table: the default. Per-atom electron delocalization (EDDB integrated over the atom) and total electron population.

Summary: three population numbers — ED (total electrons in the active subset), EDLB (localized; ED − EDDB), and EDDB (delocalized). For UHF / ROHF the α and β components are reported separately in the orbital and atomic tables.

Cyclic Delocalization Index (CDI): reported only for --pathway runs. CDI equals the smallest per-atom EDDB along the pathway — the maximum population of electrons that can delocalize uniformly around the cycle. CDI ≈ 0.9 per atom is strongly aromatic; CDI ≈ 0 is non- or anti-aromatic.

Closing line and the warning counter: the run terminates with a one-line time-stamp such as > Calculation completed on 2026-MM-DD. Total time: Xs. Whenever any non-fatal warning fired during the run, the closing line is annotated with a parenthetical count: > Calculation completed (with 1 warning!) on 2026-MM-DD. Total time: Xs. The annotation appears even in true-quiet (-q) mode where the warning text itself is suppressed — so a batch user always knows when something needs a closer look. To see the actual warning message, rerun the same command without -q. ECP-related inconsistencies (missing [Pseudo] in molden, Q-Chem GUI=2 quirk in fchk) are fatal errors — they stop the run rather than just counting; see §3 for the recovery steps.

8. Known Issues and Limitations

Tested QC programs and formats: validated against .fchk files from Gaussian 16 and Q-Chem 6 (HF / DFT and post-HF), .json and .molden files from ORCA 6 (HF / DFT and post-HF via orca_2json and orca_2mkl), .molden files from Multiwfn-converted wavefunctions, native .tmol bundles from Turbomole, and .gms output files from GAMESS-US (RHF / UHF / ROHF / MCSCF on RUNTYP = ENERGY / OPTIMIZE / SADPOINT). Cross-program agreement on benzene at cc-pVDZ is 4–5 significant figures for both MP2 and CCSD; the WF6 ECP benchmark agrees to ~3 significant figures across Gaussian / ORCA fchk / json / molden at RHF / def2-TZVPP. The built-in GFN2-xTB engine for the .xyz lane is validated against the reference tblite implementation across 81 systems covering 26 elements (H…W) including first / second / third-row transition metals, lanthanides, and main-group hypervalent compounds — final per-atom Mulliken charges agree to ~10−6e at tight SCC convergence, and total energies match Grimme’s reference xtb binary to within ~10−4 Eh.

Supported basis functions. The native lane (.fchk / .json / .molden) supports angular momenta up to k-type (l = 7), covering every standard basis set in routine use including cc-pV6Z. Overlaps are evaluated by direct Gauss-Hermite quadrature (no Obara-Saika recursion), so each shell is treated explicitly. The legacy .49 visualisation lane supports up to h-type (l = 5) only; for cc-pV6Z and similar (i / k functions) please use the matching native input format.

Post-HF density requirements. The QC program must be told to write the relaxed (or linearized) correlated density:

If the correlated block is missing, runEDDB silently falls back to the SCF reference density and the run is (correctly) classified as plain HF/DFT. The Reference type: line in the verbose-mode Stage 1 detail block reports the actual classification.

HLRL — visualisation at large basis sets. From cc-pVQZ onwards (worse at cc-pV5Z, def2-QZVPP) the rigorous NAO basis vectors pick up a small fraction of high-L (d, f, g) character from inter-orbital orthogonalization — an inherent property of the NAO algorithm, not a runEDDB artefact (NBO7 shows the same pattern). On --output-ao for native inputs, runEDDB applies a Galerkin S-metric projection on the file-side coefficients to keep the orbital pictures clean; analysis numbers are unchanged.

NBO archive (.49) — legacy lane. Population analyses (NOBD / orbital / atomic tables, summary) and all six visualisation outputs work via a parallel NBO-convention pipeline that mirrors the legacy runEDDB_nbo49 writers. The .49 path is provided for legacy compatibility only and will be removed once the native NAO driver completes its acceptance suite; for new projects, please use .fchk / .json / .molden.

Cube grid quality. The default grid (~100³ points) gives integration accuracy within ~1% for most systems. For publication images, use -c 150 or higher. The integrated density and its percentage of the analytical value are printed after each cube export for quality control.

9. Performance Tips

Parallel BOP analysis: use --ncores N to parallelize the BOP step. Speedup is nearly linear for large systems (100+ atoms).

 terminal
$ runEDDB -i large_molecule.fchk --ncores 4

GFN2-xTB on multiple cores: for the .xyz lane, --ncores N also parallelizes the GFN2-xTB SCC itself (the heaviest step for large molecules). The engine auto-sizes its thread pools to the cores you request, so requesting more cores than the job can use no longer wastes CPU. As a rule of thumb, the GFN2-xTB SCC keeps speeding up to roughly --ncores 16–32 on a few-hundred-atom system; beyond that the gains taper off (it is dominated by linear-algebra steps that do not scale indefinitely). For small molecules, --ncores 1–4 is plenty.

Analysis basis subset: the default -b nvb (Val-only) is the fastest option, since the BOP cycle scales with the number of NAOs in the working basis. -b nmb (Cor + Val) costs only marginally more (~10–30% depending on element). -b nao (Cor + Val + Ryd) is significantly slower at large basis sets — use it for completeness checks rather than routine analysis.

Cube grid sizing: the default --cube-size 100 (~1M points) is fine for routine work; for publication-quality figures and large systems we recommend at least 150. Grid evaluation can be demanding but parallelizes perfectly, so use as many cores as you can (--ncores / -n).

10. Citations

If you use runEDDB in your research, please consider to cite the following works:

  1. Electron delocalization index based on bond order orbitals.
    D.W. Szczepanik (), E.J. Zak, K. Dyduch, J. Mrozek
    Chemical Physics Letters  593 (2014) 154−159. DOI: 10.1016/j.cplett.2014.01.006.   URL 

  2. A uniform approach to the description of multicenter bonding.
    D.W. Szczepanik (), M. Andrzejak, K. Dyduch, E.J. Zak, M. Makowski, G. Mazur, J. Mrozek,
    Physical Chemistry Chemical Physics  16 (2014) 20514−20523. DOI: 10.1039/c4cp02932a.   URL 

  3. A new perspective on quantifying electron localization and delocalization in molecular systems.
    D.W. Szczepanik ()
    Computational and Theoretical Chemistry  1080 (2016) 33−37. DOI: 10.1016/j.comptc.2016.02.003.   URL 

  4. On the three-center orbital projection formalism within the electron density of delocalized bonds method.
    D.W. Szczepanik ()
    Computational and Theoretical Chemistry  1100 (2017), 13−17. DOI: 10.1016/j.comptc.2016.12.003.   URL 

  5. The electron density of delocalized bonds (EDDB) applied for quantifying aromaticity.
    D.W. Szczepanik (), M. Andrzejak, J. Dominikowska, B. Pawełek, T.M. Krygowski, H. Szatylowicz, M. Solà
    Physical Chemistry Chemical Physics  19 (2017) 28970−28981. DOI: 10.1039/c7cp06114e.   URL 

  6. Electron delocalization in planar metallacycles: Hückel or Möbius aromatic?
    D.W. Szczepanik (), M. Solà ()
    ChemistryOpen  8 (2019) 219−227. DOI: 10.1002/open.201900014.   URL 

  7. The electron density of delocalized bonds (EDDB) as a measure of local and global aromaticity
    D.W. Szczepanik (), M. Solà ()
    in
    Aromaticity: Modern Computational Methods and Applications
    Edited by I. Fernández
    Chapter 8 (pp. 259−283), Elsevier, 2021.
    ISBN: 978-0-12-822723-7
    DOI: 10.1016/B978-0-12-822723-7.00008-X.

When you analyse a plain .xyz geometry through the built-in GFN2-xTB lane, the underlying electronic-structure method is not ours — it is the work of Grimme and co-workers. Please also cite the relevant original method papers:

  1. GFN2-xTB — An Accurate and Broadly Parametrized Self-Consistent Tight-Binding Quantum Chemical Method with Multipole Electrostatics and Density-Dependent Dispersion Contributions.
    C. Bannwarth, S. Ehlert, S. Grimme
    Journal of Chemical Theory and Computation  15 (2019) 1652−1671. DOI: 10.1021/acs.jctc.8b01176.   URL 
    — the GFN2-xTB Hamiltonian used for every .xyz calculation (minimal valence vSTO-nG basis, multipole / anisotropic electrostatics, density-dependent dispersion).

  2. A generally applicable atomic-charge dependent London dispersion correction.
    E. Caldeweyher, S. Ehlert, A. Hansen, H. Neugebauer, S. Spicher, C. Bannwarth, S. Grimme
    The Journal of Chemical Physics  150 (2019) 154122. DOI: 10.1063/1.5090222.   URL 
    — the self-consistent DFT-D4 dispersion correction (charge-dependent C6 coefficients + ATM three-body term) included in the GFN2-xTB energy.

  3. High-throughput screening of spin states for transition metal complexes with spin-polarized extended tight-binding methods.
    H. Neugebauer, B. Bädorf, S. Ehlert, A. Hansen, S. Grimme
    Journal of Computational Chemistry  44 (2023) 2120−2129. DOI: 10.1002/jcc.27185.   URL 
    — the spin-polarized extension (spGFN2-xTB) used for open-shell systems (multiplicity > 1).

runEDDB’s GFN2-xTB lane implements the GFN2-xTB and DFT-D4 methods described above and is validated against Grimme’s reference implementation. Portions of it are derived from the corresponding open-source reference code — dftd4 and tblite (LGPL-3.0-or-later) and mctc-lib (Apache-2.0) — and are used and redistributed under those licenses. See the THIRD-PARTY-LICENSES.txt file included with the program for full attributions and license texts.

11. Contact

Author and main contributor: dr hab. Dariusz W. Szczepanik, prof. UJ

Contributors: dr Ouissal El Bakouri (Postdoc), mgr Mikolaj Dratwinski (PhD), mgr Pawel A. Wieczorkiewicz (PhD), mgr Joanna Zams (PhD).

Project e-mail: runeddb@gmail.com


HOME      PROJECTS      PUBLICATIONS      runEDDB      molANO      X