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]).
[1]: Tight-Binding Formalism in the Context of the PythTB Package, https://www.physics.rutgers.edu/pythtb/_downloads/915304f3240dca549efa8f491463a797/pythtb-formalism.pdf
