From c5906e7e6d08b8fdc2c3311daa3d577f9a292aae Mon Sep 17 00:00:00 2001 From: guanjihuan <34735497+guanjihuan@users.noreply.github.com> Date: Thu, 11 Nov 2021 03:36:53 +0800 Subject: [PATCH] guan-0.0.29 --- API_reference.py | 10 +-- PyPI/setup.cfg | 2 +- PyPI/src/guan/calculate_Green_functions.py | 61 ++++++++++++++++++- ...late_band_structures_and_wave_functions.py | 6 +- PyPI/src/guan/calculate_conductance.py | 43 ------------- PyPI/src/guan/calculate_scattering_matrix.py | 2 - 6 files changed, 69 insertions(+), 55 deletions(-) diff --git a/API_reference.py b/API_reference.py index 87655dc..8260a44 100644 --- a/API_reference.py +++ b/API_reference.py @@ -66,6 +66,11 @@ green_nn_n = guan.green_function_nn_n(fermi_energy, h00, h01, green_nn_n_minus, green_in_n = guan.green_function_in_n(green_in_n_minus, h01, green_nn_n) green_ni_n = guan.green_function_ni_n(green_nn_n, h01, green_ni_n_minus) green_ii_n = guan.green_function_ii_n(green_ii_n_minus, green_in_n_minus, h01, green_nn_n, green_ni_n_minus) +transfer = guan.transfer_matrix(fermi_energy, h00, h01) +right_lead_surface, left_lead_surface = guan.surface_green_function_of_lead(fermi_energy, h00, h01) +right_self_energy, left_self_energy = guan.self_energy_of_lead(fermi_energy, h00, h01) +right_self_energy, left_self_energy = self_energy_of_lead_with_h_LC_and_h_CR(fermi_energy, h00, h01, h_LC, h_CR) +green, gamma_right, gamma_left = green_function_with_leads(fermi_energy, h00, h01, h_LC, h_CR, center_hamiltonian) # calculate density of states # Source code: https://py.guanjihuan.com/calculate_density_of_states total_dos = guan.total_density_of_states(fermi_energy, hamiltonian, broadening=0.01) @@ -77,9 +82,6 @@ local_dos = guan.local_density_of_states_for_cubic_lattice_using_dyson_equation( local_dos = guan.local_density_of_states_for_square_lattice_with_self_energy_using_dyson_equation(fermi_energy, h00, h01, N2, N1, right_self_energy, left_self_energy, internal_degree=1, broadening=0.01) # calculate conductance # Source code: https://py.guanjihuan.com/calculate_conductance -transfer = guan.transfer_matrix(fermi_energy, h00, h01) -right_lead_surface, left_lead_surface = guan.surface_green_function_of_lead(fermi_energy, h00, h01) -right_self_energy, left_self_energy = guan.self_energy_of_lead(fermi_energy, h00, h01) conductance = guan.calculate_conductance(fermi_energy, h00, h01, length=100) conductance_array = guan.calculate_conductance_with_fermi_energy_array(fermi_energy_array, h00, h01, length=100) conductance = guan.calculate_conductance_with_disorder(fermi_energy, h00, h01, disorder_intensity=2.0, disorder_concentration=1.0, length=100) @@ -109,7 +111,7 @@ guan.plot(x_array, y_array, xlabel='x', ylabel='y', title='', filename='a', show guan.plot_3d_surface(x_array, y_array, matrix, xlabel='x', ylabel='y', zlabel='z', title='', filename='a', show=1, save=0, z_min=None, z_max=None) guan.plot_contour(x_array, y_array, matrix, xlabel='x', ylabel='y', title='', filename='a', show=1, save=0) -# others # Source code: https://py.guanjihuan.com/source-code/others +# others # Source code: https://py.guanjihuan.com/source-code/others guan.download_with_scihub(address=None, num=1) guan.str_to_audio(str='hello world', rate=125, voice=1, read=1, save=0, print_text=0) guan.txt_to_audio(txt_path, rate=125, voice=1, read=1, save=0, print_text=0) diff --git a/PyPI/setup.cfg b/PyPI/setup.cfg index 507ccf1..0fb0c6e 100644 --- a/PyPI/setup.cfg +++ b/PyPI/setup.cfg @@ -1,7 +1,7 @@ [metadata] # replace with your username: name = guan -version = 0.0.28 +version = 0.0.29 author = guanjihuan author_email = guanjihuan@163.com description = An open source python package diff --git a/PyPI/src/guan/calculate_Green_functions.py b/PyPI/src/guan/calculate_Green_functions.py index 89d66bd..6d14406 100644 --- a/PyPI/src/guan/calculate_Green_functions.py +++ b/PyPI/src/guan/calculate_Green_functions.py @@ -32,4 +32,63 @@ def green_function_ni_n(green_nn_n, h01, green_ni_n_minus): def green_function_ii_n(green_ii_n_minus, green_in_n_minus, h01, green_nn_n, green_ni_n_minus): green_ii_n = green_ii_n_minus+np.dot(np.dot(np.dot(np.dot(green_in_n_minus, h01), green_nn_n), h01.transpose().conj()),green_ni_n_minus) - return green_ii_n \ No newline at end of file + return green_ii_n + +def transfer_matrix(fermi_energy, h00, h01): + h01 = np.array(h01) + if np.array(h00).shape==(): + dim = 1 + else: + dim = np.array(h00).shape[0] + transfer = np.zeros((2*dim, 2*dim), dtype=complex) + transfer[0:dim, 0:dim] = np.dot(np.linalg.inv(h01), fermi_energy*np.identity(dim)-h00) + transfer[0:dim, dim:2*dim] = np.dot(-1*np.linalg.inv(h01), h01.transpose().conj()) + transfer[dim:2*dim, 0:dim] = np.identity(dim) + transfer[dim:2*dim, dim:2*dim] = 0 + return transfer + +def surface_green_function_of_lead(fermi_energy, h00, h01): + h01 = np.array(h01) + if np.array(h00).shape==(): + dim = 1 + else: + dim = np.array(h00).shape[0] + fermi_energy = fermi_energy+1e-9*1j + transfer = transfer_matrix(fermi_energy, h00, h01) + eigenvalue, eigenvector = np.linalg.eig(transfer) + ind = np.argsort(np.abs(eigenvalue)) + temp = np.zeros((2*dim, 2*dim), dtype=complex) + i0 = 0 + for ind0 in ind: + temp[:, i0] = eigenvector[:, ind0] + i0 += 1 + s1 = temp[dim:2*dim, 0:dim] + s2 = temp[0:dim, 0:dim] + s3 = temp[dim:2*dim, dim:2*dim] + s4 = temp[0:dim, dim:2*dim] + right_lead_surface = np.linalg.inv(fermi_energy*np.identity(dim)-h00-np.dot(np.dot(h01, s2), np.linalg.inv(s1))) + left_lead_surface = np.linalg.inv(fermi_energy*np.identity(dim)-h00-np.dot(np.dot(h01.transpose().conj(), s3), np.linalg.inv(s4))) + return right_lead_surface, left_lead_surface + +def self_energy_of_lead(fermi_energy, h00, h01): + h01 = np.array(h01) + right_lead_surface, left_lead_surface = surface_green_function_of_lead(fermi_energy, h00, h01) + right_self_energy = np.dot(np.dot(h01, right_lead_surface), h01.transpose().conj()) + left_self_energy = np.dot(np.dot(h01.transpose().conj(), left_lead_surface), h01) + return right_self_energy, left_self_energy + +def self_energy_of_lead_with_h_LC_and_h_CR(fermi_energy, h00, h01, h_LC, h_CR): + h_LC = np.array(h_LC) + h_CR = np.array(h_CR) + right_lead_surface, left_lead_surface = surface_green_function_of_lead(fermi_energy, h00, h01) + right_self_energy = np.dot(np.dot(h_CR, right_lead_surface), h_CR.transpose().conj()) + left_self_energy = np.dot(np.dot(h_LC.transpose().conj(), left_lead_surface), h_LC) + return right_self_energy, left_self_energy + +def green_function_with_leads(fermi_energy, h00, h01, h_LC, h_CR, center_hamiltonian): + dim = np.array(center_hamiltonian).shape[0] + right_self_energy, left_self_energy = self_energy_of_lead_with_h_LC_and_h_CR(fermi_energy, h00, h01, h_LC, h_CR) + green = np.linalg.inv(fermi_energy*np.identity(dim)-center_hamiltonian-left_self_energy-right_self_energy) + gamma_right = (right_self_energy - right_self_energy.transpose().conj())*1j + gamma_left = (left_self_energy - left_self_energy.transpose().conj())*1j + return green, gamma_right, gamma_left \ No newline at end of file diff --git a/PyPI/src/guan/calculate_band_structures_and_wave_functions.py b/PyPI/src/guan/calculate_band_structures_and_wave_functions.py index 04828eb..268d9ca 100644 --- a/PyPI/src/guan/calculate_band_structures_and_wave_functions.py +++ b/PyPI/src/guan/calculate_band_structures_and_wave_functions.py @@ -11,8 +11,7 @@ def calculate_eigenvalue(hamiltonian): if np.array(hamiltonian).shape==(): eigenvalue = np.real(hamiltonian) else: - eigenvalue, eigenvector = np.linalg.eig(hamiltonian) - eigenvalue = np.sort(np.real(eigenvalue)) + eigenvalue, eigenvector = np.linalg.eigh(hamiltonian) return eigenvalue def calculate_eigenvalue_with_one_parameter(x_array, hamiltonian_function): @@ -64,8 +63,7 @@ def calculate_eigenvalue_with_two_parameters(x_array, y_array, hamiltonian_funct ## calculate wave functions def calculate_eigenvector(hamiltonian): - eigenvalue, eigenvector = np.linalg.eig(hamiltonian) - eigenvector = eigenvector[:, np.argsort(np.real(eigenvalue))] + eigenvalue, eigenvector = np.linalg.eigh(hamiltonian) return eigenvector ## find vector with the same gauge diff --git a/PyPI/src/guan/calculate_conductance.py b/PyPI/src/guan/calculate_conductance.py index 6ce44ba..c362299 100644 --- a/PyPI/src/guan/calculate_conductance.py +++ b/PyPI/src/guan/calculate_conductance.py @@ -6,49 +6,6 @@ import numpy as np import copy from .calculate_Green_functions import * -def transfer_matrix(fermi_energy, h00, h01): - h01 = np.array(h01) - if np.array(h00).shape==(): - dim = 1 - else: - dim = np.array(h00).shape[0] - transfer = np.zeros((2*dim, 2*dim), dtype=complex) - transfer[0:dim, 0:dim] = np.dot(np.linalg.inv(h01), fermi_energy*np.identity(dim)-h00) - transfer[0:dim, dim:2*dim] = np.dot(-1*np.linalg.inv(h01), h01.transpose().conj()) - transfer[dim:2*dim, 0:dim] = np.identity(dim) - transfer[dim:2*dim, dim:2*dim] = 0 - return transfer - -def surface_green_function_of_lead(fermi_energy, h00, h01): - h01 = np.array(h01) - if np.array(h00).shape==(): - dim = 1 - else: - dim = np.array(h00).shape[0] - fermi_energy = fermi_energy+1e-9*1j - transfer = transfer_matrix(fermi_energy, h00, h01) - eigenvalue, eigenvector = np.linalg.eig(transfer) - ind = np.argsort(np.abs(eigenvalue)) - temp = np.zeros((2*dim, 2*dim), dtype=complex) - i0 = 0 - for ind0 in ind: - temp[:, i0] = eigenvector[:, ind0] - i0 += 1 - s1 = temp[dim:2*dim, 0:dim] - s2 = temp[0:dim, 0:dim] - s3 = temp[dim:2*dim, dim:2*dim] - s4 = temp[0:dim, dim:2*dim] - right_lead_surface = np.linalg.inv(fermi_energy*np.identity(dim)-h00-np.dot(np.dot(h01, s2), np.linalg.inv(s1))) - left_lead_surface = np.linalg.inv(fermi_energy*np.identity(dim)-h00-np.dot(np.dot(h01.transpose().conj(), s3), np.linalg.inv(s4))) - return right_lead_surface, left_lead_surface - -def self_energy_of_lead(fermi_energy, h00, h01): - h01 = np.array(h01) - right_lead_surface, left_lead_surface = surface_green_function_of_lead(fermi_energy, h00, h01) - right_self_energy = np.dot(np.dot(h01, right_lead_surface), h01.transpose().conj()) - left_self_energy = np.dot(np.dot(h01.transpose().conj(), left_lead_surface), h01) - return right_self_energy, left_self_energy - def calculate_conductance(fermi_energy, h00, h01, length=100): right_self_energy, left_self_energy = self_energy_of_lead(fermi_energy, h00, h01) for ix in range(length): diff --git a/PyPI/src/guan/calculate_scattering_matrix.py b/PyPI/src/guan/calculate_scattering_matrix.py index e695a0c..6cc4d5d 100644 --- a/PyPI/src/guan/calculate_scattering_matrix.py +++ b/PyPI/src/guan/calculate_scattering_matrix.py @@ -5,8 +5,6 @@ import numpy as np import copy from .calculate_Green_functions import * -from .calculate_conductance import * - def if_active_channel(k_of_channel): if np.abs(np.imag(k_of_channel))<1e-6: