xy_coordinates Module

Defining real space 2D coordinates for honeycomb lattice sites


Uses

  • module~~xy_coordinates~~UsesGraph module~xy_coordinates xy_coordinates assert_m assert_m module~xy_coordinates->assert_m module~hex_coordinates hex_coordinates module~xy_coordinates->module~hex_coordinates module~hex_layout hex_layout module~xy_coordinates->module~hex_layout module~hex_coordinates->assert_m

Used by

  • module~~xy_coordinates~~UsedByGraph module~xy_coordinates xy_coordinates module~honeyplots honeyplots module~honeyplots->module~xy_coordinates module~xy_neighbors xy_neighbors module~honeyplots->module~xy_neighbors module~honeytools honeytools module~honeytools->module~xy_coordinates module~honeytools->module~xy_neighbors module~xy_neighbors->module~xy_coordinates

Variables

Type Visibility Attributes Name Initial
integer, private, parameter :: N = 6
real(kind=8), private, parameter :: PI = 4d0*atan(1d0)

Interfaces

public interface xy_print

  • private impure elemental subroutine xy_class_print(S, unit, quiet)

    Pretty print of xy coordinates in static arrays

    Arguments

    Type IntentOptional Attributes Name
    class(xy), intent(in) :: S
    integer, intent(in), optional :: unit

    default = $stdout

    logical, intent(in), optional :: quiet

    default = .false.

  • private impure elemental subroutine xy_lattice_print(S, unit, quiet)

    Pretty print of xy coordinates in dynamic arrays

    Arguments

    Type IntentOptional Attributes Name
    type(xy_lattice), intent(in) :: S
    integer, intent(in), optional :: unit

    default = $stdout

    logical, intent(in), optional :: quiet

    default = .false.


Derived Types

type, public ::  xy

Base type for 2D points

Components

Type Visibility Attributes Name Initial
real(kind=8), public :: x
real(kind=8), public :: y

Type-Bound Procedures

generic, public :: operator(==) => eq_xy
generic, public :: operator(/=) => neq_xy
procedure, private :: eq_xy
procedure, private :: neq_xy

type, public, extends(xy) ::  xy_site

A 2D point extension for inequivalent sites

Components

Type Visibility Attributes Name Initial
real(kind=8), public :: x
real(kind=8), public :: y
character(len=1), public :: label

A or B (sublattice)

integer, public :: key = 0

for lattice lookup

Type-Bound Procedures

generic, public :: operator(==) => eq_xy
generic, public :: operator(/=) => neq_xy

type, public ::  xy_lattice

Wrapper type for storing dynamically sized collections of lattice sites

Components

Type Visibility Attributes Name Initial
type(xy_site), public, allocatable :: site(:)
integer, public :: size

Type-Bound Procedures

generic, public :: operator(==) => eq_lattice
generic, public :: operator(/=) => neq_lattice
procedure, public :: push_back
procedure, private :: eq_lattice
procedure, private :: neq_lattice

Functions

public pure elemental function xy_distance(A, B) result(d)

Polymorphic euclidean distance for xy class

Arguments

Type IntentOptional Attributes Name
class(xy), intent(in) :: A
class(xy), intent(in) :: B

Return Value real(kind=8)

public pure elemental function hex2site(layout, H, label) result(site)

Convert hex coordinates to real 2D lattice sites [returning the xy coordinates for the unique unit-cell, sites "A" and "B": this of course does not account for border effects, so would be suitable only for simple, periodized, systems, which for now are out of scope.] Actual real-space layout has to be specified by passing a (scalar) unit_cell object.

Arguments

Type IntentOptional Attributes Name
type(unit_cell), intent(in) :: layout
type(hex), intent(in) :: H
character(len=1), intent(in) :: label

A or B

Return Value type(xy_site)

public pure elemental function hex2corner(layout, H) result(corner)

Convert hex coordinates to real 2D space [returning the xy coordinates for the hexagon corners, as appropriatiely wrapped in the "xy_lattice" type] Actual real-space layout has to be specified by passing a (scalar) unit_cell object.

Arguments

Type IntentOptional Attributes Name
type(unit_cell), intent(in) :: layout
type(hex), intent(in) :: H

Return Value type(xy_lattice)

public pure elemental function hex2center(layout, H) result(center)

Convert hex coordinates to real 2D space [returning the xy coordinates of the hexagon center] Actual real-space layout has to be specified by passing a (scalar) unit_cell object.

Arguments

Type IntentOptional Attributes Name
type(unit_cell), intent(in) :: layout
type(hex), intent(in) :: H

Return Value type(xy)

public pure function hex2lattice(layout, hexagons) result(lattice)

Generate a type(xy_lattice) object from any given hex array, provided a suitable layout (unit-cell)

Arguments

Type IntentOptional Attributes Name
type(unit_cell), intent(in) :: layout
type(hex), intent(in) :: hexagons(:)

Return Value type(xy_lattice)

public pure elemental function get_sublattice(lattice, label) result(sublattice)

Extract sublattice, given a lattice and a label ("A" or "B")

Arguments

Type IntentOptional Attributes Name
type(xy_lattice), intent(in) :: lattice
character(len=1), intent(in) :: label

A or B

Return Value type(xy_lattice)

public pure function xy_ordered_union(A, B) result(C)

An ordered-keys set union for xy_lattices. It compares the lattice sites, looking at their x and y coordinates only: equal x,y entries are inserted in C just once, so to make it a set. A and B must be sets too. The sublattice labels are preserved in the process, assuming that two sites with the same x,y pertain to the same sublattice. For A the keys are assumed to be 1:size(A), and preserved as such (with an assertion). The keys of B are instead discarded, so to be replaced by size(A):size(C), which thus amounts to have result C uniquely indexed as 1, 2, ... , size(A) + size(B). This allows building consistently indexed matrices to act on the lattice array, such as real-space tight-binding hamiltonians. The keys would then be used to index other real space quantities, such as LDOS, Chern marker, local magnetization, and so on.

Arguments

Type IntentOptional Attributes Name
type(xy_lattice), intent(in) :: A
type(xy_lattice), intent(in) :: B

Return Value type(xy_lattice)

private pure function ith_corner_offset(layout, i) result(offset)

Compute the i-th offset vector connecting the hexagon center to its corners, returning a (scalar) xy coordinate

Arguments

Type IntentOptional Attributes Name
type(unit_cell), intent(in) :: layout
integer, intent(in) :: i

Return Value type(xy)

private pure elemental function eq_xy(A, B) result(isequal)

tolerance equality for xy class [absolute tol hard-coded to 1d-12] -> probably to improve...

Arguments

Type IntentOptional Attributes Name
class(xy), intent(in) :: A
class(xy), intent(in) :: B

Return Value logical

private pure elemental function neq_xy(A, B) result(notequal)

tolerance inequality for xy class [absolute tol hard-coded to 1d-12] -> probably to improve...

Arguments

Type IntentOptional Attributes Name
class(xy), intent(in) :: A
class(xy), intent(in) :: B

Return Value logical

private pure elemental function eq_lattice(A, B) result(isequal)

polymorphic equality overload for xy_lattice type

Arguments

Type IntentOptional Attributes Name
class(xy_lattice), intent(in) :: A
class(xy_lattice), intent(in) :: B

Return Value logical

private pure elemental function neq_lattice(A, B) result(notequal)

polymorphic inequality overload for xy_lattice type

Arguments

Type IntentOptional Attributes Name
class(xy_lattice), intent(in) :: A
class(xy_lattice), intent(in) :: B

Return Value logical

private pure elemental function xy_norm(A) result(r)

polymorphic euclidean norm for xy class

Arguments

Type IntentOptional Attributes Name
class(xy), intent(in) :: A

Return Value real(kind=8)


Subroutines

private pure subroutine push_back(vec, val)

Poor man implementation of a dynamic array, a là std::vector (but without preallocation and smart doubling...)

Arguments

Type IntentOptional Attributes Name
class(xy_lattice), intent(inout) :: vec
type(xy_site), intent(in) :: val

private impure elemental subroutine xy_class_print(S, unit, quiet)

Pretty print of xy coordinates in static arrays

Arguments

Type IntentOptional Attributes Name
class(xy), intent(in) :: S
integer, intent(in), optional :: unit

default = $stdout

logical, intent(in), optional :: quiet

default = .false.

private impure elemental subroutine xy_lattice_print(S, unit, quiet)

Pretty print of xy coordinates in dynamic arrays

Arguments

Type IntentOptional Attributes Name
type(xy_lattice), intent(in) :: S
integer, intent(in), optional :: unit

default = $stdout

logical, intent(in), optional :: quiet

default = .false.