jdftx icon indicating copy to clipboard operation
jdftx copied to clipboard

Eliashberg spectral function α2F in terms of the mode-resolved coupling strengths λqν

Open dijasila opened this issue 1 year ago • 3 comments

Hello everyone, Is it possible in the current implementation to calculate the The Eliashberg spectral function α2F d in terms of the mode-resolved coupling strengths λqν and the phonon frequencies? is there a simple way of doing that? All I want is to precisely know which peak in the α2F correspond to which specific phonon mode..

Best Regards

dijasila avatar Jan 15 '24 03:01 dijasila

Please see the script within http://jdftx.org/EphMatrixElements.html that computes the Eliashberg spectral function on teh route to compute resistivity. In particular see this section of code from the penultimate block:

# Term combining matrix elements and velocity factors, summed over bands:
    w1 = np.exp(-0.5*((E1-mu)/Esigma)**2) * (1./(Esigma*np.sqrt(2*np.pi)))
    w2 = np.exp(-0.5*((E2-mu)/Esigma)**2) * (1./(Esigma*np.sqrt(2*np.pi)))
    term = np.einsum('ka,Kb,kKab,kKxab->kKx', w1, w2, vFactor, np.abs(g)**2)
    # Histogram by phonon frequency:
    F += np.histogram(omegaPh.flatten(), omegaPhBins, weights=term.flatten())[0]

In that term is calculated with indices kKx where x is phonon mode, and this is flattened and summed over in the next line. You can easily modify this section to separately histogram contributions from each phonon mode to get the mode-resolved spectral function.

Best, Shankar

shankar1729 avatar Jan 15 '24 19:01 shankar1729

Hi Shankar, Apologies for my poor understanding of the code that calculates Eliashberg functions. Can you please elaborate a bit more on your response? So what you are saying is I should process the 'term' seperately for each phonon mode and not take a sum as is in line "F += np.histogram(omegaPh.flatten(), omegaPhBins, weights=term.flatten())[0]"?

something like this.... term = np.einsum('ka,Kb,kKab,kKxab->kKx', w1, w2, vFactor, np.abs(g)**2) # Histogram by phonon frequency: F = np.histogram(omegaPh.flatten(), omegaPhBins, weights=term.flatten())[0] F[iblock] += F

dijasila avatar Jan 16 '24 00:01 dijasila

No, something like this (omitting the unmodified block of code in the middle, showing with ...):

F = np.zeros((nModes, len(omegaPhMid)))  # spectral function for each phonon mode
print('Sampling transport spectral function:', end=' ')
for iBlock in range(Nblocks):
    # Select blocks of k-points from above set:

    ...

    # Histogram by phonon frequency, separately for each phonon mode:
    for iMode in range(nModes):
        F[iMode] += np.histogram(omegaPh[..., iMode].flatten(), omegaPhBins, weights=term[..., iMode].flatten())[0]
    print(iBlock+1, end=' '); sys.stdout.flush()

where the additional dimension of F is each phonon mode. You have to of course modify the subsequent plotting and analysis routines accordingly. (I have not tested this explicitly; just sketching out what would need to be done.)

shankar1729 avatar Jan 16 '24 21:01 shankar1729