# What is the best program to manipulate numerical DFT wavefunctions to calculate custom matrix elements?

For pedagogical reasons, I am looking for ways to calculate quantum-mechanical quantities such as $$\langle m | \dot{m}\rangle, \langle m | \dot{n}\rangle, \langle m | \ddot{n}\rangle$$ using wavefunctions $$m,n$$ output from DFT calculations (Quantum Espresso, specifically). The overhead dots denote derivatives. Ideally, I would be able to calculate the Berry phase for custom loops as well.

So far, I have tried using some Python tools (z2Pack, PythTB) to try and post-process wannier90 output, but as a beginner, it doesn't seem as if these tools allow us to calculate our own matrix elements. They seem to be just for intra-level quantities such as $$\langle m | \dot{m}\rangle$$.

Does anyone more-experienced have any advice on what existing software might be easiest to work with to calculate these custom-defined quantities? I would rather try messing with some software that is more likely to support these kinds of calculations. PythTB seems to be one of the best options I have, but Python may not be the best for large datasets (compared to Fortran, etc). Thank you for your time.

Edit: The issue also seems to be that several post-processing software tend to use the numerical method where one takes $$\arg$$ of some product of complex phases corresponding to each $$k$$ point. However to do something like $$\langle m | \dot{n}\rangle$$, it might be better to use a central difference method to carry out the derivative of the wavefunction (as opposed to an established discretized method for the Berry phase, as in section 4.5 in Ref [1]).

You can use QE.6xxx with the support of the hdf5 library. To realize that purpose, you should add the following command when you compile QE:

--with-hdf5=yes

or take a look at the official guide.

Then the saved wavefunction can be manipulated with the python package h5py.

import h5py