From ff550b33c35e53416d20284db2e542db237b6999 Mon Sep 17 00:00:00 2001 From: Florian Knoop Date: Fri, 19 Jan 2024 18:16:01 +0100 Subject: [PATCH] canonical_configuration | add flag to displace imaginary modes --- src/canonical_configuration/main.f90 | 4 +-- src/canonical_configuration/options.f90 | 6 ++++ .../type_forceconstant_secondorder.f90 | 3 +- .../type_forceconstant_secondorder_aux.f90 | 29 +++++++++++++++---- 4 files changed, 33 insertions(+), 9 deletions(-) diff --git a/src/canonical_configuration/main.f90 b/src/canonical_configuration/main.f90 index c5aaca52..5eb97468 100644 --- a/src/canonical_configuration/main.f90 +++ b/src/canonical_configuration/main.f90 @@ -188,9 +188,9 @@ program canonical_configuration else invert = .false. end if - call fcss%initialize_cell(p, uc, fc, opts%temperature, opts%zpm, .false., opts%mindist, mw, imode=imode, invert=invert) + call fcss%initialize_cell(p, uc, fc, opts%temperature, opts%zpm, .false., opts%mindist, mw, imode=imode, invert=invert, imaginary=opts%imaginary) else - call fcss%initialize_cell(p, uc, fc, opts%temperature, opts%zpm, .false., opts%mindist, mw) + call fcss%initialize_cell(p, uc, fc, opts%temperature, opts%zpm, .false., opts%mindist, mw, imaginary=opts%imaginary) end if ! dump to file diff --git a/src/canonical_configuration/options.f90 b/src/canonical_configuration/options.f90 index db68e36c..c2606354 100644 --- a/src/canonical_configuration/options.f90 +++ b/src/canonical_configuration/options.f90 @@ -26,6 +26,7 @@ module options real(r8) :: dielcutoff2 = -lo_huge real(r8) :: dielcutoff3 = -lo_huge logical :: modes = .false. + logical :: imaginary = .false. contains procedure :: parse end type @@ -118,6 +119,10 @@ subroutine parse(opts) call cli%add(switch='--modes', & help='Print displacements for every individual mode.', hidden=.true., & required=.false., act='store_true', def='.false.', error=lo_status) + call cli%add(switch='--imaginary', hidden=.true., & + help='Displace imaginary modes as if they were positive', & + required=.false., act='store_true', def='.false.', error=lo_status) + if (lo_status .ne. 0) stop cli_manpage @@ -153,6 +158,7 @@ subroutine parse(opts) call cli%get(switch='-dc2', val=opts%dielcutoff2) call cli%get(switch='-dc3', val=opts%dielcutoff3) call cli%get(switch='--modes', val=opts%modes) + call cli%get(switch='--imaginary', val=opts%imaginary) ! Convert input to atomic units right away opts%maximum_frequency = opts%maximum_frequency*lo_frequency_THz_to_hartree diff --git a/src/libolle/type_forceconstant_secondorder.f90 b/src/libolle/type_forceconstant_secondorder.f90 index 4ed7c70a..a501c7f3 100644 --- a/src/libolle/type_forceconstant_secondorder.f90 +++ b/src/libolle/type_forceconstant_secondorder.f90 @@ -190,7 +190,7 @@ module subroutine fake_forceconstant(fc, uc, ss, debye_temperature, maximum_freq real(r8), intent(in), optional :: maximum_frequency integer, intent(in), optional :: verbosity end subroutine - module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, closest_distance, mw, nosync, imode, invert) + module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, closest_distance, mw, nosync, imode, invert, imaginary) class(lo_forceconstant_secondorder), intent(inout) :: fcss type(lo_crystalstructure), intent(inout) :: ss type(lo_crystalstructure), intent(in) :: uc @@ -203,6 +203,7 @@ module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, logical, intent(in), optional :: nosync integer, intent(in), optional :: imode logical, intent(in), optional :: invert + logical, intent(in), optional :: imaginary end subroutine module subroutine setsumtozero(fc) class(lo_forceconstant_secondorder), intent(inout) :: fc diff --git a/src/libolle/type_forceconstant_secondorder_aux.f90 b/src/libolle/type_forceconstant_secondorder_aux.f90 index d56d22ff..9750f5b0 100644 --- a/src/libolle/type_forceconstant_secondorder_aux.f90 +++ b/src/libolle/type_forceconstant_secondorder_aux.f90 @@ -363,7 +363,7 @@ subroutine set_values(alpha, fc) end subroutine !> use the harmonic model to initialize a cell -module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, closest_distance, mw, nosync, imode, invert) +module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, closest_distance, mw, nosync, imode, invert, imaginary) !> force constant for this (super) cell class(lo_forceconstant_secondorder), intent(inout) :: fcss !> supercell to be thermally populated @@ -388,7 +388,10 @@ module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, integer, intent(in), optional :: imode !> use negative mode amplitude so that u_i = \sum_s -Q_s X_si) logical, intent(in), optional :: invert + !> displace imaginary modes instead of discarding + logical, intent(in), optional :: imaginary real(r8) :: inv_prefactor + logical :: imaginary_ ! Not sure about save attribute here. type(lo_mersennetwister), save :: tw @@ -402,6 +405,12 @@ module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, end if end if + if (present(imaginary)) then + imaginary_ = imaginary + else + imaginary_ = .false. + end if + init: block ! Seed rng if needed if (tw%initialized .eqv. .false.) then @@ -435,13 +444,21 @@ module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, ! Set the amplitudes setamplitude: block integer :: i + real(r8) :: omega_tmp + do i = 1, fcss%na*3 - if (fcss%omega(i) .gt. lo_freqtol) then + if (imaginary_) then + ! imaginary modes are supposed to be displaced instead of frozen + omega_tmp = abs(fcss%omega(i)) + else + omega_tmp = fcss%omega(i) + end if + if (omega_tmp .gt. lo_freqtol) then ! Choose quantum or classical statistics if (quantum) then - fcss%amplitudes(i) = sqrt((2*lo_planck(temperature, fcss%omega(i)) + 1)*0.5_r8/fcss%omega(i)) + fcss%amplitudes(i) = sqrt((2*lo_planck(temperature, omega_tmp) + 1)*0.5_r8/omega_tmp) else - fcss%amplitudes(i) = sqrt(lo_kb_hartree*temperature)/fcss%omega(i) + fcss%amplitudes(i) = sqrt(lo_kb_hartree*temperature)/omega_tmp end if else ! set to zero for acoustic modes @@ -474,7 +491,7 @@ module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, do j = 1, 3 l = l + 1 ss%u(j, a1) = +fcss%amplitudes(imode)*ss%invsqrtmass(a1)*fcss%eigenvectors(l, imode) - ss%v(j, a1) = -fcss%amplitudes(imode)*ss%invsqrtmass(a1)*fcss%eigenvectors(l, imode)*fcss%omega(imode) + ss%v(j, a1) = -fcss%amplitudes(imode)*ss%invsqrtmass(a1)*fcss%eigenvectors(l, imode)*abs(fcss%omega(imode)) end do end do else @@ -487,7 +504,7 @@ module subroutine initialize_cell(fcss, ss, uc, fc, temperature, quantum, exact, do j = 1, 3 l = l + 1 ss%u(j, a1) = ss%u(j, a1) + fcss%amplitudes(i)*ss%invsqrtmass(a1)*x1*fcss%eigenvectors(l, i) - ss%v(j, a1) = ss%v(j, a1) - fcss%amplitudes(i)*ss%invsqrtmass(a1)*x2*fcss%eigenvectors(l, i)*fcss%omega(i) + ss%v(j, a1) = ss%v(j, a1) - fcss%amplitudes(i)*ss%invsqrtmass(a1)*x2*fcss%eigenvectors(l, i)*abs(fcss%omega(i)) end do end do end do modeloop