WTP provides some IO capabilities for interfacing with Quantum Espresso and Wannier90.

This documentation represents my best effort on figuring out the file formats. It does not serve as a reference, but it's probably as good as any. If there is a difference from the documentation of Quantum Espresso/Wannier90, you should rely on their documentation instead (if you can find them :P).

`wfcX.dat`

files

The `wfcX.dat`

files are used by Quantum Espresso to store the wave functions in binary format. It seems that they are moving to HDF5, but I couldn't find any documentation on that just yet as of 2021.

### File format

`WTP.WFC`

β TypeThis structure provides an interface to the .wfc files.

The format is documented in: format of wfc files.

The relevant part is copied here for convenience:

`INTEGER :: ik`

k-point index (1 to number of k-points)`REAL(8) :: xk(3)`

k-point coordinates`INTEGER :: ispin`

spin index for LSDA case: ispin=1 for spin-up, ispin=2 for spin-down for unpolarized or non-colinear cases, ispin=1 always`LOGICAL :: gamma_only`

if .true. write or read only half of the plane waves`REAL(8) :: scalef`

scale factor applied to wavefunctions`INTEGER :: ngw`

number of plane waves (PW)`INTEGER :: igwx`

max number of PW (may be larger than ngw, not sure why)`INTEGER :: npol`

number of spin states for PWs: 2 for non-colinear case, 1 otherwise`INTEGER :: nbnd`

number of wavefunctions`REAL(8) :: b1(3), b2(3), b3(3)`

primitive reciprocal lattice vectors`INTEGER :: mill(3,igwx)`

miller indices:`h=mill(1,i), k=mill(2,i), l=mill(3,i)`

the i-th PW has wave vector`(k+G)(:)=xk(:)+h*b1(:)+k*b2(:)+ l*b3(:)`

`COMPLEX(8) :: evc(npol*igwx,nbnd)`

wave functions in the PW basis set The first index runs on PW components, the second index runs on band states. For non-colinear case, each PW has a spin component first igwx components have PW with up spin, second igwx components have PW with down spin

Some variables have been renamed in the package to exorcise the acronyms.

`ik`

->`i_kpoint`

`xk`

->`k_coordinates`

`ispin`

->`i_spin`

`scalef`

->`scale_factor`

`ngw`

->`n_planewaves`

`igwx`

->`max_n_planewaves`

`npol`

->`n_polerizations`

`nbnd`

->`n_band`

`b1, b2, b3`

->`b_1, b_2, b_3`

`WTP.WFC`

β Method`WFC(wave_function_filename)`

Load the metadata of a wave function from a wfc.dat (terrible filenames) file. This does not load the wave functions.

This is not guaranteed to work when Quantum Espresso is compiled with a Fortran compiler that is not `gfort`

. Fortran does not enforce a standard on the size of it's primitive types in binary files.

Example:

Assume that you are in `docs`

directory (change the path otherwise).

```
julia> path_to_si = "../test/test_data/test_5";
julia> wave_function = WFC(joinpath(path_to_si, "si.save/wfc1.dat"));
julia> wave_function.n_band, wave_function.gamma_only, wave_function.n_planewaves
(4, false, 537)
```

### Read a `.save`

Directory

`WTP.wave_functions_from_directory`

β Function`wave_functions_from_directory(save_dir)`

Read all the wave function files in an directory for metadata. This does not load the wave functions.

The `wfcX.dat`

files only gives you the content of a single wave function, but what you really want is probably reading the entire `.save`

directory and convert them to orbitals on grids.

```
julia> wave_functions_list = wave_functions_from_directory(joinpath(path_to_si, "si.save"));
julia> length(wave_functions_list)
64
```

`WTP.orbital_set_from_save`

β Function```
orbital_set_from_save(
wave_functions_list,
[domain_scaling_factor=2,]
)
```

Load the wavefunctions and parse them into orbitals.

The `domain_scaling_factor`

is possibly either from Nyquist-Shannon sampling theorem. or from the energy cutoff for the density. Settting this to 2 gives the same FFT as QE.

```
julia> u = orbital_set_from_save(wave_functions_list);
julia> length(elements(u))
64
```

`WTP.single_orbital_from_wave_functions`

β Functionsingle*orbital*from*wave*functions(w, l, k, n)

Use this if you want to load a single orbital. Mostly likely you will want `orbital_set_from_save`

to load them all at once.

Create an orbital on the nth band from a WFC object w. The orbital will be defined on a reciprocal lattice l.

k is the kpoint, which can be out of the first Brillouin zone. If k is outside the first Brillouin zone, we will fold the kpoint into the first Brillouin zone and translate the wave function in on the reciprocal lattice (phase shift in real lattice).

The orbital is represented as values on the reciprocal lattice.

`WTP.load_evc!`

β Function`load_evc!(wfc, [bands=1:wfc.n_band])`

Load the evc (frequency space wave-functions) into wave_functions. You probably shouldn't have to touch this.

### K-point map

`WTP.i_kpoint_map`

β Function`i_kpoint_map(wave_function_list)`

Construct a mapping between the mystical i_kpoints and kpoints from a list of WFC objects. Also figure out the Brillouin zone. For technical reasons, these two things are difficult to do separately.

This mapping may be different from what is used for .mmn and .amn files.

```
julia> k_map, brillouin_zone = i_kpoint_map(wave_functions_list);
julia> k_map[1]
GridVector{BrillouinZone3D}:
coefficients: [0, 0, 0]
```

`UNKXXXXX.1`

files

`WTP.UNK`

β TypeThe `UNKXXXXX.1`

files are produced by `pw2wannier90.x`

to show orbitals in real space.

This structure provides an interface to the UNK files.

`nx, ny, nz`

: Number of grid points in x, y, z direction.`k`

: The kpoint at which the unk orbital is defined.`n_band`

: The number of bands stored.`psi_r`

: Each row of this matrix corresponds to a band. Each column corresponds to a grid point (not sure which one yet).`filename`

: The name of the file from which the UNK object is constructed.

`WTP.UNK`

β Method`UNK(unk_filename)`

Load a unk file named unk_filename as a UNK object.

The wave functions are read such that each row is an orbital to avoid transposing in memory for SCDM.

This takes 10s to load a 1.4G unk file using 12 threads and an NVMe drive. Since the unk files are too clumsy to cinlude in VCS, we will not provide the files.

Example:

`wave_function = UNK(joinpath(path_to_si, "unk/UNK00001.1"))`

`WTP.wannier_from_unk_dir`

β Function```
wannier_from_unk_dir(
unk_dir,
wave_functions_list,
[domain_scaling_factor=2,]
)
```

Load the wavefunctions and parse them into orbitals. Since the `UNKXXXXX.1`

does not carry much information, we need the `wave_functions_list`

from the `.save`

directory as well.

Example:

`u_real = wannier_from_unk_dir(joinpath(path_to_si, "unk"), wave_functions_list)`

`WTP.single_orbital_from_unk`

β Function`single_orbital_from_unk(unk, h, k, n)`

Use this if you want to load a single orbital. Mostly likely you will want `wannier_from_unk_dir`

to load them all at once.

Create an orbital on the nth band from a UNK structure unk. The orbital will be defined on a homecell h.

k is the kpoint, which can be out of the first Brillouin zone.

*The kpoint is not folded and the wave function is not phase shifted*. Currently, this is done at a later stage after the fft. TODO: Phase shift the orbital.

The orbital is represented as values on the homecell.

`.mmn`

Files

`WTP.MMN`

β TypeThis structure provides an interface to the `MMN`

files.

- The first line should tell a joke. The parser won't parse unless tickled.

The format of the `MMN`

files can be found in the user guide of `wannier90.x`

.

Ξ-point:

`wannier90.x`

uses MMN files also for the X, Y, Z matrices of gamma point calculations. This is nowhere documented.

The $g$-vectors (taken from wannier90 user guide):

The last three integers specify the $G$ vector, in reciprocal lattice units, that brings the k-point specified by the second integer, and that thus lives inside the first Brillouin zone, to the actual $k + b$ we need.

`WTP.MMN`

β Method`MMN(mmn_filename, [ad_hoc_gamma = false])`

Load a MMN file. Set `ad_hoc_gamma = true`

if it is a gamma-only calculation so that the gamma-trick shenanigans are dealt with.

```
julia> mmn = MMN(joinpath(path_to_si, "output/pw2wan/si.mmn"))
julia> mmn.n_kpoint, mmn.n_neighbor
(64, 8)
```

<!β

Missing docstring for `NeighborIntegral(::MMN, ::Dict{Int64, KPoint})`

. Check Documenter's build log for details.

β>

For Ξ-point calculations, the neighbor integrals are between Ξ-point of the first Brillouin zone and that of other (most likely adjacent) Brillouin zones.

`.amn`

Files

`WTP.AMN`

β TypeThe `.amn`

files are used by Wannier90 for storing the gauge transform.

`WTP.AMN`

β Method`AMN(amn_filename)`

Create an AMN object from a file.

```
julia> amn = AMN(joinpath(path_to_si, "output/pw2wan/si.amn"));
julia> amn.n_kpoint
64
```

`WTP.Gauge`

β Method`Gauge(grid, amn, k_map, [orthonormalization = true])`

Create a gauge from an `AMN`

object and a `k_map`

.

An `AMN`

object also just stores the raw data, and we need to convert this to a `Gauge`

to use it. Like `.mmn`

files, we need a mapping of k-points as well as the Brillouin zone to construct the gauge, which can then be indexed with k-points.

```
julia> U = Gauge(brillouin_zone, amn, k_map);
julia> U[brillouin_zone[0, 0, 0]][:, :]
4Γ4 Matrix{ComplexF64}:
-0.424021+0.0687693im -0.530967+0.0861146im -0.469845+0.0762009im -0.540275+0.0876234im
0.00775601-0.282049im -0.3875-0.232102im 0.714604+0.092279im -0.246713+0.369214im
0.837758+0.116721im -0.437746+0.134205im -0.142503-0.20453im -0.103362-0.045631im
0.0694477+0.124822im 0.516605+0.173543im 0.152318-0.411003im -0.694672+0.0889086im
```