Output (HDF5)

HDF5 is the recommended file format for writing and storing DCA++'s output. As it is a binary format, it cannot be viewed in a normal text editor, however, the h5py package provides a simple python interface to HDF5 files.

show_hdf5.py

The python script show_hdf5.py provides a first overview of the content stored in an HDF5 file. Running it on a typical output file from a main_dca run, called dca_sp.hdf5 in the example, produces the following list ([...] indicates omissions):

$ python show_hdf5.py dca_sp.hdf5

CT-AUX-SOLVER-functions
[...]
DCA-loop-functions
[...]
domains
domains/CLUSTER
domains/CLUSTER/MOMENTUM_SPACE
domains/CLUSTER/MOMENTUM_SPACE/basis
domains/CLUSTER/MOMENTUM_SPACE/basis/data
domains/CLUSTER/MOMENTUM_SPACE/basis/equal-size
domains/CLUSTER/MOMENTUM_SPACE/basis/size
domains/CLUSTER/MOMENTUM_SPACE/elements
domains/CLUSTER/MOMENTUM_SPACE/elements/data
domains/CLUSTER/MOMENTUM_SPACE/elements/equal-size
domains/CLUSTER/MOMENTUM_SPACE/elements/size
domains/CLUSTER/MOMENTUM_SPACE/super-basis
domains/CLUSTER/MOMENTUM_SPACE/super-basis/data
domains/CLUSTER/MOMENTUM_SPACE/super-basis/equal-size
domains/CLUSTER/MOMENTUM_SPACE/super-basis/size
domains/CLUSTER/REAL_SPACE
[...]
domains/frequency-domain
domains/frequency-domain/elements
domains/frequency-domain/indices
[...]
functions
functions/Self_Energy
functions/Self_Energy/data
functions/Self_Energy/domain-sizes
functions/Self_Energy/name
functions/Self_Energy-stddev
functions/Self_Energy-stddev/data
functions/Self_Energy-stddev/domain-sizes
functions/Self_Energy-stddev/name
[...]
parameters
parameters/CT-AUX
parameters/CT-AUX/additional-time-measurements
parameters/CT-AUX/expansion-parameter-K
parameters/CT-AUX/initial-configuration-size
parameters/CT-AUX/initial-matrix-size
parameters/CT-AUX/max-submatrix-size
parameters/CT-AUX/neglect-Bennett-updates
parameters/DCA
[...]

Reading and plotting a dca::func::function

The most relevant data sets in a DCA++ output file correspond to objects stored as dca::func::functions and the domains they are defined on. These data sets can easily be loaded into NumPy arrays in order to further analyze and plot the functions.

Example

This C++ code snippet defines the cluster self-energy Σ_c (see dca_data.hpp):

using namespace dca::phys;

using b = func::dmn_0<domains::electron_band_domain>;  // 1 element for a simple square lattice
using s = func::dmn_0<domains::electron_spin_domain>;  // 2 elements (up, down)
using wn = func::dmn_0<domains::frequency_domain>;     // Fermionic Matsubara frequency

// Cluster momenta
using K = func::dmn_0<domains::cluster_domain<double, 2, domains::CLUSTER,
	                                          domains::MOMENTUM_SPACE, domains::BRILLOUIN_ZONE>>;

// Cluster self-energy
func::function<std::complex<double>, func::dmn_variadic<b, s, b, s, K, wn>> self_energy("Self_Energy");

Assuming the self-energy, defined as above, has been written to a file called dca_sp.hdf5 together with its subdomains, k_DCA and wn, we can plot it as a function of Matsubara frequency using the python code below.
Note that only applying the slice operator [start:stop:step] actually loads the data set into a NumPy array, and [:] loads the whole data set. Also notice the shape of the self-energy array. The order of the subdomains is reversed, while real and imaginary part correspond to the last dimension.

import h5py
from matplotlib import pyplot as plt

# The HDF5 file containing the data.
filename = 'dca_sp.hdf5'

# Open the file.
data = h5py.File(filename,'r')

# Matsubara frequency mesh
wn = data['domains/frequency-domain/elements'] [:]

# DCA cluster momenta (not needed for the plot)
K = data['domains/CLUSTER/MOMENTUM_SPACE/elements/data'][:]

# Cluster self-energy
# shape: [wn, K, spin2, band2, spin1, band1, real/imag]
self_energy = data['functions/Self_Energy/data']

# Plot the cluster self-energy vs. Matsubara frequency for the first cluster momentum ([0.,0.]).
# All band and spin indices are set to zero.
K_ind = 0
b1 = b2 = s1 = s2 = 0

plt.plot(wn, self_energy[:, K_ind, s2, b2, s1, b1, 0], label='real')
plt.plot(wn, self_energy[:, K_ind, s2, b2, s1, b1, 1], label='imag')
plt.legend(loc='best')
plt.show()