From 33891747c22fa296d0258f6e9439cc34dd0d849a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alo=C3=AFs=20Castellano?= Date: Thu, 8 May 2025 09:38:23 +0200 Subject: [PATCH 1/4] Add a modecoupling flag to remove the real part of the self energy at omega=0 --- src/lineshape/options.f90 | 12 ++++++++++++ src/lineshape/phonondamping_generation.f90 | 8 ++++++++ 2 files changed, 20 insertions(+) diff --git a/src/lineshape/options.f90 b/src/lineshape/options.f90 index 16ae51b5..cd8c76fe 100644 --- a/src/lineshape/options.f90 +++ b/src/lineshape/options.f90 @@ -20,6 +20,7 @@ module options logical :: isotopescattering = .false. logical :: thirdorder = .false. logical :: fourthorder = .false. + logical :: mct = .false. logical :: slightsmearing = .false. integer :: integrationtype = -lo_hugeint logical :: readiso = .false. @@ -129,6 +130,10 @@ subroutine parse(opts) help='Consider four-phonon contributions to the real part of the self-energy.', hidden=.true., & required=.false., act='store_true', def='.false.', error=lo_status) if (lo_status .ne. 0) stop + call cli%add(switch='--modecoupling', switch_ab='-mct', & + help='Remove the static contribution to the self-energy in the mode-coupling theory.', & + required=.false., act='store_true', def='.false.', error=lo_status) + if (lo_status .ne. 0) stop call cli%add(switch='--nondiagonal', & help='Consider non-diagonal contributions to the self-energy.', hidden=.true., & required=.false., act='store_true', def='.false.', error=lo_status) @@ -233,6 +238,7 @@ subroutine parse(opts) call cli%get(switch='--no_thirdorder_scattering', val=dumlog) opts%thirdorder = .not. dumlog call cli%get(switch='--fourthorder', val=opts%fourthorder) + call cli%get(switch='--modecoupling', val=opts%mct) call cli%get(switch='--nondiagonal', val=dumlog) opts%diagonal = .not. dumlog call cli%get(switch='--minsmear', val=opts%minsmear) @@ -329,6 +335,12 @@ subroutine parse(opts) end if end if + ! If we are in the mode-coupling approach, the real part four phonon makes no sense + if (opts%mct .and. opts%fourthorder) then + write(*, *) 'There is no real part from the fourth order in the mode-coupling approximation' + stop + end if + ! Not really options opts%timereversal = .true. opts%slightsmearing = .false. diff --git a/src/lineshape/phonondamping_generation.f90 b/src/lineshape/phonondamping_generation.f90 index bb9110f7..6d2e737c 100644 --- a/src/lineshape/phonondamping_generation.f90 +++ b/src/lineshape/phonondamping_generation.f90 @@ -247,8 +247,16 @@ subroutine generate(se, qpoint, qdir, wp, uc, fc, fct, fcf, ise, isf, qp, dr, op call mem%deallocate(xs, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) call mem%deallocate(z0, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) call mem%deallocate(y0, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + + ! If we are in the mode-coupling approximation, we need to remove the static contribution + if (opts%mct) then + do imode=1, se%n_mode + se%re_3ph(:, imode) = se%re_3ph(:, imode) - se%re_3ph(1, imode) + end do + end if end block kktransform call tmr%tock('Kramers-Kronig transformation') + end if ! Normalize the spectral functions From 9177301e3ee80cbdcc6f61006ab598597fa62a8f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alo=C3=AFs=20Castellano?= Date: Thu, 8 May 2025 09:39:38 +0200 Subject: [PATCH 2/4] redid a stopping message --- src/lineshape/options.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/lineshape/options.f90 b/src/lineshape/options.f90 index cd8c76fe..7be1a128 100644 --- a/src/lineshape/options.f90 +++ b/src/lineshape/options.f90 @@ -337,7 +337,7 @@ subroutine parse(opts) ! If we are in the mode-coupling approach, the real part four phonon makes no sense if (opts%mct .and. opts%fourthorder) then - write(*, *) 'There is no real part from the fourth order in the mode-coupling approximation' + write(*, *) 'There is no real part from the fourth order in the mode-coupling theory' stop end if From 0c1527043ea56e3c927134568ebfaf62ee7567b5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alo=C3=AFs=20Castellano?= Date: Tue, 13 May 2025 15:14:47 +0200 Subject: [PATCH 3/4] adding orthogonal time correction to lineshape --- src/lineshape/phonondamping_generation.f90 | 59 ++++++++++++++++++++++ 1 file changed, 59 insertions(+) diff --git a/src/lineshape/phonondamping_generation.f90 b/src/lineshape/phonondamping_generation.f90 index 6d2e737c..870fd9f4 100644 --- a/src/lineshape/phonondamping_generation.f90 +++ b/src/lineshape/phonondamping_generation.f90 @@ -256,7 +256,66 @@ subroutine generate(se, qpoint, qdir, wp, uc, fc, fct, fcf, ise, isf, qp, dr, op end if end block kktransform call tmr%tock('Kramers-Kronig transformation') + end if + + ! Correct the memory matrix to get back to the orthogonal time one + if (opts%mct .and. se%thirdorder_scattering) then + kktransform2: block + real(r8) :: pref = 2.0_r8/lo_pi + complex(r8), dimension(:), allocatable :: z0 + complex(r8) :: eta + real(r8), dimension(:), allocatable :: x, xs, y0 + real(r8) :: xp, dlx + integer :: imode, ie, ctr + + call mem%allocate(x, se%n_energy, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + call mem%allocate(xs, se%n_energy, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + call mem%allocate(z0, se%n_energy, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + call mem%allocate(y0, se%n_energy, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + x = se%energy_axis ! copy of x-axis + xs = x**2 ! x^2, precalculated + dlx = x(2) - x(1) ! prefactor for riemann integral + eta = lo_imag*(x(2) - x(1))*1E-8_r8 ! small imaginary thing to avoid divergence + + do imode=1, se%n_mode + se%im_3ph(2:, imode) = (x(2:)**2 * se%im_3ph(2:, imode)) / ((x(2:) - se%re_3ph(2:, imode))**2 + se%im_3ph(2:, imode)**2) + end do + + ctr = 0 + se%re_3ph = 0.0_r8 + do imode = 1, se%n_mode + do ie = 1, se%n_energy + ctr = ctr + 1 + if (mod(ctr, mw%n) .ne. mw%r) cycle + xp = x(ie) + y0 = se%im_3ph(:, imode)*x + ! To anyone reading this code: + ! The correct way to compute the principal value would be + ! z0 = (xp + eta)**2 - xs + ! However, it doesn't matter for the Im -> Re transform because we + ! don't really hit 0. + z0 = xp**2 - xs + eta + y0 = real(y0/z0, r8) + y0(1) = y0(1)*0.5_r8 + y0(se%n_energy) = y0(se%n_energy)*0.5_r8 + se%re_3ph(ie, imode) = sum(y0)*dlx + end do + end do + call mw%allreduce('sum', se%re_3ph) + se%re_3ph = se%re_3ph*pref + + call mem%deallocate(x, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + call mem%deallocate(xs, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + call mem%deallocate(z0, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + call mem%deallocate(y0, persistent=.false., scalable=.false., file=__FILE__, line=__LINE__) + + ! If we are in the mode-coupling approximation, we need to remove the static contribution + do imode=1, se%n_mode + se%re_3ph(:, imode) = se%re_3ph(:, imode) - se%re_3ph(1, imode) + end do + end block kktransform2 + call tmr%tock('Kramers-Kronig transformation') end if ! Normalize the spectral functions From 654a2acb4ce154e941339fda35099102483d2ca9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alo=C3=AFs=20Castellano?= Date: Fri, 16 Jan 2026 14:44:53 +0100 Subject: [PATCH 4/4] Gives a better name for the option --- docs/program/lineshape.md | 4 ++++ src/lineshape/options.f90 | 6 +++--- 2 files changed, 7 insertions(+), 3 deletions(-) diff --git a/docs/program/lineshape.md b/docs/program/lineshape.md index ac5884b9..3480b262 100644 --- a/docs/program/lineshape.md +++ b/docs/program/lineshape.md @@ -86,6 +86,10 @@ Optional switches: default value .false. Read the q-point mesh from file. To generate a q-mesh file, see the genkpoints utility. +* `--remove_static_selfenergy` + default value .false. + Remove the static contribution to the self-energy as in the mode-coupling theory. + * `--help`, `-h` Print this help message diff --git a/src/lineshape/options.f90 b/src/lineshape/options.f90 index 7be1a128..2bd48b32 100644 --- a/src/lineshape/options.f90 +++ b/src/lineshape/options.f90 @@ -130,8 +130,8 @@ subroutine parse(opts) help='Consider four-phonon contributions to the real part of the self-energy.', hidden=.true., & required=.false., act='store_true', def='.false.', error=lo_status) if (lo_status .ne. 0) stop - call cli%add(switch='--modecoupling', switch_ab='-mct', & - help='Remove the static contribution to the self-energy in the mode-coupling theory.', & + call cli%add(switch='--remove_static_selfenergy', & + help='Remove the static contribution to the self-energy as in the mode-coupling theory.', & required=.false., act='store_true', def='.false.', error=lo_status) if (lo_status .ne. 0) stop call cli%add(switch='--nondiagonal', & @@ -238,7 +238,7 @@ subroutine parse(opts) call cli%get(switch='--no_thirdorder_scattering', val=dumlog) opts%thirdorder = .not. dumlog call cli%get(switch='--fourthorder', val=opts%fourthorder) - call cli%get(switch='--modecoupling', val=opts%mct) + call cli%get(switch='--remove_static_selfenergy', val=opts%mct) call cli%get(switch='--nondiagonal', val=dumlog) opts%diagonal = .not. dumlog call cli%get(switch='--minsmear', val=opts%minsmear)