smaal corrections

README.md

@@ -55,7 +55,7 @@ The following methods are defined within the class:
 
 These routines are based on the publication "Computation of the surface electron-energy-loss spectrum in specular geometry for an arbitrary plane-stratified medium" by P. Lambin, J.-P. Vigneron, and A. A. Lucas, in the Journal "Computer Physics Communications 60, 351-64(1990)".
 
-The code is modified to comply with Fortran90 and wraped to python function. See the example in calcHREELS1.py for NiO(001). The parameters for epsilon_infinity "eps": 5.25, as well as frequency and width for the TO phonon,  "wTO": [393.7], "gTO": [10.8], and the LO phonon, "wLO": [584.7], "gLO": [10.8], need to be specified.
+The code is modified to comply with Fortran90 and wraped to python functions. See the example in Examples/calcHREELS1.py for NiO(001). The parameters for epsilon_infinity "eps": 5.25, as well as frequency and width for the TO phonon,  "wTO": [393.7], "gTO": [10.8], and the LO phonon, "wLO": [584.7], "gLO": [10.8], need to be specified.
 
 Complex calculations for perovskite oxides are provided in the examples calcHREELS2.py and calcHREELS3.py.
 

libhreels/calcHREELS20.py

@@ -73,9 +73,6 @@ def addDrude(wLOPlasma, gLOPlasma, material, gTOPlasma='None'):
     if gTOPlasma == 'none':
         gTOPlasma = gLOPlasma
     newMaterial = deepcopy(material)
-    # To select the Drude response, the Lambin Fortran code requires a negative first oscillator frequency
-    # which is then interpreted as plasma frequency.
-    # The values of the former LO parameter are irrelevant (but need to be provided). 
     try:
         if len(newMaterial['wTO']) > 0:
             newMaterial['wTO'] += [0.]
@@ -198,52 +195,19 @@ def myMain():
     import os
 
     x = np.linspace(-100.,1000,400)
-    x2 = np.linspace(-100.,1000,2400)
-    
-
     material = {'eps': 4.,
-                'wTO': [0],  'gTO': [20], 'wLO': [598.7], 'gLO': [20],
-                }
-    
-    print("Ag imported as: ", importMaterials('Ag'))
+                'wTO': [0],  'gTO': [20], 'wLO': [598.7], 'gLO': [20]}
 
     film1 = lambin(film=[[material,10000.]])
     film1.temperature = 100
-    #xs, spectrum = film1.calcSurfaceLoss(x)
-    xs, spectrum = film1.calcHREELS(x,normalized=False,areanormalized=False)
-    xs1, spectrum1 = film1.calcHREELS(x2,normalized=False,areanormalized=False)
-    xs2, spectrum2 = film1.calcHREELS(x,areanormalized=True)
-    xs3, spectrum3 = film1.calcHREELS(x2,areanormalized=True)
-
-    norm_test=integrate.simps(spectrum, dx=xs[2]-xs[1])/(max(xs)-min(xs))
-    norm_test1=integrate.simps(spectrum1, dx=xs1[2]-xs1[1])/(max(xs1)-min(xs1))
-    norm_test2=integrate.simps(spectrum2, dx=xs2[2]-xs2[1])/(max(xs2)-min(xs2))
-    norm_test3=integrate.simps(spectrum3, dx=xs3[2]-xs3[1])/(max(xs3)-min(xs3))
-    print(norm_test,norm_test1,norm_test2,norm_test3)
-    print(xs[2]-xs[1],xs1[2]-xs1[1],xs2[2]-xs2[1],xs3[2]-xs3[1])
-    
-    plt.plot(xs[:-1],spectrum[:-1], label='normalized=Flase, '+str(len(x))+' points, E0 at '+str(max(spectrum)))
-    plt.plot(xs1[:-2],spectrum1[:-2], label='normalized=Flase, '+str(len(x2))+' points, E0 at '+str(max(spectrum1)))
-    plt.plot(xs2[:-10],spectrum2[:-10],label='area normalized=True, '+str(len(x))+' points, E0 at '+str(max(spectrum2)))
-    plt.plot(xs3[:-20],spectrum3[:-20],label='area normalized=True, '+str(len(x2))+' points, E0 at '+str(max(spectrum3)))
-
-    print('spec2/spec0=', max(spectrum2)/max(spectrum))
-    print('spec3/spec0=', max(spectrum3)/max(spectrum))
-    print('spec0/spec2=', max(spectrum)/max(spectrum2))
-    print('spec0/spec3=', max(spectrum)/max(spectrum3))
-
-    # pureDrude = {'eps': 4.,
-    #             'wTO': [-400],  'gTO': [-20], 'wLO': [0], 'gLO': [0],
-    #             }
-    # film2 = lambin(film=[[pureDrude,10000.]])
-    # film2.temperature = 100
-    # xs, spectrum = film2.calcSurfaceLoss(x)
-    # plt.plot(xs,spectrum,'g-.', label='pure Drude 400')
-
+    xs, spectrum = film1.calcHREELS(x,normalized=True,areanormalized=False)
 
+    plt.plot(xs[:-1],spectrum[:-1], label='normalized=False')
     plt.ylabel('Surface Loss')
     plt.xlabel('Energy Loss (cm$^{-1}$)')
-    plt.legend(title=r'Plasma frequency')
+    plt.ylim(bottom=0)
+    plt.title('Bulk plasmon at 600 cm-1')
+    # plt.legend(title=r'Plasma frequency')
 
     plt.text(0.99, 0.01,os.path.basename(__file__), fontsize=10, ha='right', va='bottom', transform=plt.gcf().transFigure)
     output_filename = os.path.splitext(__file__)[0] + '.png'

libhreels/dielectrics20.py

 # Testing dielectric models and their surface loss function
 import matplotlib.pyplot as plt
 import numpy as np
-from numpy import real, imag, sqrt
 from scipy import constants
 
-
-def doping2chargecarrierdensity(doping): #use doping in percentage by mass, not volume
-    '''Returns the charge carrier density for a given doping in percentage of mass.'''
-    return float(doping*3.31664067935517E+020*1E6) # doping[%] * charge carrier/cm^3 * 1E6 = n_e [m^-3] 
-
-def doping2plasmaFrequency(doping,epsInfinity=1.): #use doping in percentage by mass, not volume
-    '''Returns the plasma frequency for a given doping in percentage of mass.'''
+def plasmaFrequency(n): 
+    '''Returns the plasma frequency for a given doping.'''
+    # Note that the function parameter is actually n/m, the charge carrier density 
+    # over the band mass in units of electron mass in vacuum.
     eps_0 = constants.value("vacuum electric permittivity")
     e = constants.value('elementary charge')
     m_e = constants.value('electron mass')
-    c = constants.value('speed of light in vacuum') 
-    eps_r = epsInfinity
-    n_e = doping2chargecarrierdensity(doping)
-
-    w_P = np.sqrt(n_e*e*e/(eps_r*eps_0*m_e)) #Hz
+    w_P = np.sqrt(n*e*e/(eps_0*m_e)) #Hz
     w_P = w_P /1E+12 * 33.35641 #Hz -> THz -> cm^-1
     return w_P 
     
-
-def doping2surfacePlasma(doping,epsInfinity=1.):
-    '''Returns the surface plasma frequency in cm-1 for a doping 
-       given by the argument in cm-3'''
-    return np.sqrt(doping)*1.8817885780819758e-06/np.sqrt(1+epsInfinity)
+# def doping2chargecarrierdensity(doping): #use doping in percentage by mass, not volume
+#     '''Returns the charge carrier density for a given doping in percentage.'''
+#     return float(doping*3.31664067935517E+020*1E6) # doping[%] * charge carrier/cm^3 * 1E6 = n_e [m^-3] 
+
+# def doping2plasmaFrequency(doping): 
+#     '''Returns the plasma frequency for a given doping.'''
+#     n = doping2chargecarrierdensity(doping)
+#     return plasmaFrequency(n) 
+    
 
 def loss(eps):
     '''Returns the loss function for a given eps'''
-    return np.imag(eps) #HHE changed former sign
+    return np.imag(-1/eps) 
 
 def surfaceLoss(eps):
     '''Returns the surface loss function for a given eps'''
@@ -56,7 +52,7 @@ def plotDielectrics(x,eps, title=" ", plot_show=True):
     axs[0].legend()
     axs[0].set_ylabel('Dielectric Function')
 
-    axs[1].plot(x, surfaceLoss(eps), lineStyle='-', label='Loss function')    
+    axs[1].plot(x, surfaceLoss(eps), linestyle='-', label='Loss function')    
     axs[1].set_ylabel('Surface Loss Function')
     axs[1].set_ylim([-0.10,None])
 
@@ -97,7 +93,7 @@ class drude:
         self.gamma_p = gamma_p
         self.gamma_0 = gamma_0
 
-        print("called drude with w_P, g_P, g_0:", omega_p, gamma_p, gamma_0)
+        #print("called drude with w_P, g_P, g_0:", omega_p, gamma_p, gamma_0)
         
 
     def eps(self,w):
@@ -126,7 +122,7 @@ class cole_davidson:
         self.gamma_0 = gamma_0
         self.beta = beta
 
-        print("called coledavidson with w_P, g_P, g_0, beta:", omega_p, gamma_p, gamma_0, beta)
+        #print("called coledavidson with w_P, g_P, g_0, beta:", omega_p, gamma_p, gamma_0, beta)
     
 
     def eps(self,w):
@@ -161,7 +157,7 @@ class cole_cole:
         self.gamma_0 = gamma_0
         self.beta = beta
 
-        print("called colecole with w_P, g_P, g_0, beta:", omega_p, gamma_p, gamma_0, beta)
+        #print("called colecole with w_P, g_P, g_0, beta:", omega_p, gamma_p, gamma_0, beta)
 
     def eps(self,w):
         newW = np.where(w==0, 0.000567894, w)   # Avoid zeros to avoid devision by zero 

pyproject.toml

 [tool.poetry]
 name = "libhreels"
-version = "2.0.0"
+version = "2.0.1"
 description = "Handling, simulating, and plotting HREELS and Auger spectroscopy data"
 authors = ["Wolf Widdra <wolf.widdra@physik.uni-halle.de>"]
 include = ["*./libhreels/*"]