Skip to content

Commit

Permalink
Pre-commit fixes
Browse files Browse the repository at this point in the history
  • Loading branch information
PhilippRue committed May 10, 2024
1 parent 0118b47 commit e591719
Show file tree
Hide file tree
Showing 4 changed files with 248 additions and 255 deletions.
148 changes: 76 additions & 72 deletions aiida_kkr/tools/tools_STM_scan.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@
##############################################################################
# combine impurty clusters


def convert_to_imp_cls(host_structure, imp_info):
"""
convert imp info to rcls form
Expand Down Expand Up @@ -230,77 +231,76 @@ def create_combined_potential_node_cf(add_position, host_remote, imp_potential_n


#####################################################################
# STM pathfinder
# STM pathfinder


def STM_pathfinder(host_structure):
from aiida_kkr.tools import find_parent_structure
from ase.spacegroup import Spacegroup

"""This function is used to help visualize the scanned positions
"""This function is used to help visualize the scanned positions
and the symmetries that are present in the system """

"""
inputs:
inputs:
host_struture : RemoteData : The Remote data contains all the information needed to create the path to scan
outputs:
outputs:
struc_info : Dict : Dictionary containing the structural information of the film
matrices : Array : Array containing the matrices that generate the symmetries of the system
"""

def info_creation(structure):
from ase.spacegroup import get_spacegroup
# List of the Bravais vectors
# List of the Bravais vectors
vec_list = structure.cell.tolist()

# Find the Bravais vectors that are in plane vectors
plane_vectors = {'plane_vectors':[], 'space_group':''}
plane_vectors = {'plane_vectors': [], 'space_group': ''}
for vec in vec_list:
# Is this sufficient to find all the in-plane vectors?
if vec[2] == 0:
plane_vectors['plane_vectors'].append(vec[:2])

space_symmetry = get_spacegroup(ase_struc)
plane_vectors['space_group'] = space_symmetry.no

return plane_vectors

def symmetry_finder(struc_info):
# Here we get the symmetry operations that are possible
symmetry_matrices = Spacegroup(struc_info['space_group'])

# Reduce the dimensionality, we only want the 2D matrices
matrices = []
for element in symmetry_matrices.get_rotations():
matrices.append(element[:2, :2])

# Uniqueness of the elements must be preserved
unique_matrices = []
for matrix in matrices:
if not any(np.array_equal(matrix, m) for m in unique_matrices):
unique_matrices.append(matrix)

return unique_matrices

struc = find_parent_structure(host_structure)
# clone the structure since it has already been saved in AiiDA and cannot be modified
supp_struc = struc.clone()

# If the structure is not periodic in every direction we force it to be.
supp_struc.pbc = (True, True, True)
# ASE struc

# ASE struc
ase_struc = supp_struc.get_ase()

# Structural informations are stored here
struc_info = info_creation(ase_struc)

# The structural informations are then used to find the symmetries of the system
symm_matrices = symmetry_finder(struc_info)

return struc_info, symm_matrices

return struc_info, symm_matrices


@engine.calcfunction
def STM_pathfinder_cf(host_structure):
Expand All @@ -318,50 +318,52 @@ def STM_pathfinder_cf(host_structure):


def lattice_generation(x_len, y_len, rot, vec):
import math

import math
"""
x_len : int : value to create points between - x and x.
y_len : int : value to create points between - y and y.
rot : list : list of the rotation matrices given by the symmetry of the system.
vec : list : list containing the two Bravais vectors.
"""

# Here we create a grid in made of points whic are the linear combination of the lattice vectors
lattice_points = []
for i in range(-x_len, x_len+1):

for i in range(-x_len, x_len + 1):
lattice_points_col = []
for j in range(-y_len, y_len+1):
p = [i*x + j*y for x,y in zip(vec[0], vec[1])]
for j in range(-y_len, y_len + 1):
p = [i * x + j * y for x, y in zip(vec[0], vec[1])]
lattice_points_col.append(p)
lattice_points.append(lattice_points_col)

# Eliminiatio of the symmetrical sites
points_to_eliminate = []
for i in range(-x_len, x_len+1):
for j in range(-y_len, y_len+1):

for i in range(-x_len, x_len + 1):
for j in range(-y_len, y_len + 1):
if lattice_points[i][j][0] >= 0 and lattice_points[i][j][1] >= 0:
for element in rot[1:]:
point = np.dot(element, lattice_points[i][j])
if point[0] >= 0 and point[1] >=0:
if point[0] >= 0 and point[1] >= 0:
continue
else:
points_to_eliminate.append(point)

points_to_scan = []

for i in range(-x_len, x_len + 1):
for j in range(-y_len, y_len + 1):
eliminate = False
for elem in points_to_eliminate:
# Since there can be some numerical error in the dot product we use the isclose function
if (math.isclose(elem[0], lattice_points[i][j][0], abs_tol=1e-4) and math.isclose(elem[1], lattice_points[i][j][1], abs_tol=1e-4)):
if (
math.isclose(elem[0], lattice_points[i][j][0], abs_tol=1e-4) and
math.isclose(elem[1], lattice_points[i][j][1], abs_tol=1e-4)
):
eliminate = True
if not eliminate:
points_to_scan.append(lattice_points[i][j])

return points_to_eliminate, points_to_scan


Expand All @@ -373,48 +375,50 @@ def lattice_plot(plane_vectors, symm_vec, symm_matrices, grid_length_x, grid_len
from aiida_kkr.tools import lattice_generation
import matplotlib.pyplot as plt
from matplotlib.lines import Line2D
origin = np.array([[0, 0],[0, 0]])

origin = np.array([[0, 0], [0, 0]])
# Generation of the points to plot
unused , used = lattice_generation(grid_length_x, grid_length_y, symm_matrices, plane_vectors)
unused, used = lattice_generation(grid_length_x, grid_length_y, symm_matrices, plane_vectors)

# Plotting of the points
for element in unused:
plt.scatter(element[0], element[1], marker='s', s=130, c='#33638DFF')

for element in used:
plt.scatter(element[0], element[1], marker='D', s=130, c='#FDE725FF')
# Plot of the crystal symmetry directions, tag must be activated.
plt.scatter(element[0], element[1], marker='D', s=130, c='#FDE725FF')

# Plot of the crystal symmetry directions, tag must be activated.
if symm_vec:
import numpy.linalg as lin

for element in symm_matrices:
eig_val, eig_vec = lin.eig(element)

for element in eig_vec:
plt.quiver(*origin, element[0], element[1], alpha=1, color='#B8DE29FF', angles='xy', scale_units='xy', scale=1)

plt.quiver(
*origin, element[0], element[1], alpha=1, color='#B8DE29FF', angles='xy', scale_units='xy', scale=1
)

# Plot of the Bravais lattice
for element in plane_vectors:
plt.quiver(*origin, element[0], element[1], color='#3CBB75FF', angles='xy', scale_units='xy', scale=1)

legend_elements = [
Line2D([0], [0], color='#33638DFF', lw=2, label='Unscanned Sites', marker='s'),
Line2D([0], [0], color='#FDE725FF', lw=2, label='Scanned Sites', marker='D'),
Line2D([0], [0], color='#3CBB75FF', lw=2, label='Bravais lattice'),
Line2D([0], [0], color='#33638DFF', lw=2, label='Unscanned Sites', marker='s'),
Line2D([0], [0], color='#FDE725FF', lw=2, label='Scanned Sites', marker='D'),
Line2D([0], [0], color='#3CBB75FF', lw=2, label='Bravais lattice'),
]
plt.legend(handles=legend_elements, bbox_to_anchor=(0.75, -0.15))

plt.title('Lattice plot and symmetry directions')
plt.ylabel('y direction')
plt.xlabel('x direction')
#plt.xticks(np.arange(-grid_length, grid_length, float(plane_vectors[0][0])))
#plt.set_cmap(cmap)
plt.grid(linestyle='--')
plt.grid(linestyle='--')
plt.show()


##########################################################
# find linear combination coefficients

Expand All @@ -423,31 +427,31 @@ def find_linear_combination_coefficients(plane_vectors, vectors):
from operator import itemgetter
"""This helper function takes the planar vectors and a list of vectors
and return the coefficients in the base of the planar vectors"""

# Formulate the system of equations Ax = b
A = np.vstack((plane_vectors[0], plane_vectors[1])).T

# Solve the system of equations using least squares method
data = []
for element in vectors:
b = element
b = element
# We use the least square mean error procedure to evaulate the units of da and db
# lstsq returns: coeff, residue, rank, and singular value
# We only need the coefficient.
data.append(np.linalg.lstsq(A, b, rcond=None))

indices = []
for element in data:
supp = []
for elem in element[0]:
# Here we round to an integer, this is because of the numerical error
# which is present inside the calculation.
# Here we round to an integer, this is because of the numerical error
# which is present inside the calculation.
supp.append(round(elem))
indices.append(supp)

# Before returning the indices, we reorder them first from the lowest to the highest valued
# on the x axis and then from the lowest to the highest on the y axis.
indices = sorted(indices, key=itemgetter(0,1))
# on the x axis and then from the lowest to the highest on the y axis.

indices = sorted(indices, key=itemgetter(0, 1))

return indices
Loading

0 comments on commit e591719

Please sign in to comment.