From bbf2e8745ac7918bee04afd120998c16a54a4466 Mon Sep 17 00:00:00 2001
From: Karl-Michael Schindler <karl-michael.schindler@web.de>
Date: Thu, 29 Aug 2024 11:12:02 +0200
Subject: [PATCH] really fix of the Lorentz term.
---
source/WfW seteps/myEels20-seteps.f90 | 15 ++++++++++-----
tests/seteps/setepsdriver-WFW.f90 | 2 ++
2 files changed, 12 insertions(+), 5 deletions(-)
diff --git a/source/WfW seteps/myEels20-seteps.f90 b/source/WfW seteps/myEels20-seteps.f90
index 495f9ff..44d3ef6 100755
--- a/source/WfW seteps/myEels20-seteps.f90
+++ b/source/WfW seteps/myEels20-seteps.f90
@@ -28,7 +28,12 @@ subroutine seteps(nos, osc_size, osc, epsinf, wn, nLayer, eps)
j = 0
do l = 1, nLayer ! loop over different thin film layers
+ ! The numbering does not work with a mix of Lorentz and Kurosawa/extended Drude terms
+ ! Either Lorentz with no Drude or Kurosawa and extended Drude terms
m = nos(l)/2 ! m number of TO modes = offset to reach the LO mode in the joint TO-LO list
+ if (osc(1,1) < 0.0) then
+ m = nos(l)
+ endif
nomi = dcmplx(1.0d0, 0.0d0)
deno = dcmplx(1.0d0, 0.0d0)
addeps = dcmplx(0.0d0, 0.0d0)
@@ -37,21 +42,21 @@ subroutine seteps(nos, osc_size, osc, epsinf, wn, nLayer, eps)
do k = 1, m ! loop over all TO modes
j = j + 1
- if (osc(1,j) > 0.) then ! positive TO mode: 'Kurosawa' form: _Multiplicative_ phonon mode
+ if (osc(1,j) > 0.0) then ! positive TO mode: 'Kurosawa' form: _Multiplicative_ phonon mode
b = wn/osc(1, j+m)
nomi = nomi * osc(1, j+m)**2 * (1.0 - b * dcmplx(b, osc(3,j+m)/osc(1, j+m)) )
deno =deno * (osc(1,j)**2 - wn * dcmplx( wn, osc(3,j) ) )
-
- else if (osc(1,j) < 0.) then! Negative TO mode means: _Additive_ Lorentz oscillator with Q
+ j = j + m ! we have already looped over the LO modes, therefore increase the index
+ else if (osc(1,j) < 0.0) then! Negative TO mode means: _Additive_ Lorentz oscillator with Q
gamma = abs(osc(1,j)) * osc(3,j) ! osc(3,j) is lambda. gamma = omega_TO * lambda
addlorentz = addlorentz + osc(1,j)**2 * osc(2,j) /dcmplx(osc(1,j)**2 - wn2, -wn * gamma) ! Sign of imaginary part changed (WFW)
- else ! osc(1,j) = 0 -> it is a Drude term
+ else ! osc(1,j) = 0 -> it is a double damping extended Drude term
addeps = addeps - dcmplx(osc(1,j+m)**2, wn*(osc(3,j)-osc(3,j+m))) /dcmplx(wn2, wn*osc(3,j))
+ j = j + m ! we have already looped over the LO modes, therefore increase the index
end if
enddo
- j = j+m ! we have already looped over the LO modes, therefore increase the index
eps(l) = epsinf(l) * (nomi / deno + addeps) + addlorentz ! brackets changed by HHe 230915
enddo
debugFirstRun = .false.
diff --git a/tests/seteps/setepsdriver-WFW.f90 b/tests/seteps/setepsdriver-WFW.f90
index e0bf0a9..32a3de7 100644
--- a/tests/seteps/setepsdriver-WFW.f90
+++ b/tests/seteps/setepsdriver-WFW.f90
@@ -58,6 +58,8 @@ program setepsdriver
endif
jos = jos + 1
read(11, *) (osc(k, jos), k = 1, 3)
+! make the first parameter negative as a sign this is a Lorentz oscillator
+ osc(1,jos) = -osc(1,jos)
enddo
endif
enddo
--
GitLab