init_sites_tilted Subroutine

private subroutine init_sites_tilted(this)

Type Bound

tilted

Arguments

Type IntentOptional Attributes Name
class(tilted) :: this

Contents

Source Code


Source Code

    subroutine init_sites_tilted(this)
        class(tilted) :: this
        character(*), parameter :: this_routine = "init_sites_tilted"

        integer :: temp_array(-this%length(1):this%length(2), &
                              -this%length(1):this%length(2) + 1)
        integer :: up(-this%length(1):this%length(2), &
                      -this%length(1):this%length(2) + 1)
        integer :: down(-this%length(1):this%length(2), &
                        -this%length(1):this%length(2) + 1)
        integer :: left(-this%length(1):this%length(2), &
                        -this%length(1):this%length(2) + 1)
        integer :: right(-this%length(1):this%length(2), &
                         -this%length(1):this%length(2) + 1)
        integer :: right_ul(-this%length(1):this%length(2), &
                            -this%length(1):this%length(2) + 1)
        integer :: right_ur(-this%length(1):this%length(2), &
                            -this%length(1):this%length(2) + 1)
        integer :: right_dl(-this%length(1):this%length(2), &
                            -this%length(1):this%length(2) + 1)
        integer :: right_dr(-this%length(1):this%length(2), &
                            -this%length(1):this%length(2) + 1)
        integer :: right_rr(-this%length(1):this%length(2), &
                            -this%length(1):this%length(2) + 1)
        integer :: right_ll(-this%length(1):this%length(2), &
                            -this%length(1):this%length(2) + 1)
        integer :: up_ul(-this%length(1):this%length(2), &
                         -this%length(1):this%length(2) + 1)
        integer :: up_ur(-this%length(1):this%length(2), &
                         -this%length(1):this%length(2) + 1)
        integer :: up_dl(-this%length(1):this%length(2), &
                         -this%length(1):this%length(2) + 1)
        integer :: up_dr(-this%length(1):this%length(2), &
                         -this%length(1):this%length(2) + 1)
        integer :: up_rr(-this%length(1):this%length(2), &
                         -this%length(1):this%length(2) + 1)
        integer :: up_ll(-this%length(1):this%length(2), &
                         -this%length(1):this%length(2) + 1)
        integer :: down_ul(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: down_ur(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: down_dl(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: down_dr(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: down_rr(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: down_ll(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: left_ul(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: left_ur(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: left_dl(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: left_dr(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: left_rr(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: left_ll(-this%length(1):this%length(2), &
                           -this%length(1):this%length(2) + 1)
        integer :: i, j, k, l, pbc, temp_neigh(4), k_min, k_max, offset, k_vec(3), m
        integer :: right_nn, left_nn, up_nn, down_nn, pbc_1(2), pbc_2(2), r_vec(3)
        integer, allocatable :: neigh(:)
        integer, allocatable :: order(:)
        ! convention of lattice storage:
        !
        !   2 5
        ! 1 3 6 8
        !   4 7
        ! update: we also want to have non-square tilted clusters.. how would
        ! that work. eg a 10-site 2x3 cluster:
        !   2 5
        ! 1 3 6 9
        !   4 7 10
        !     8
        ! can i also do a 1x2 4-site tilted? like:
        ! 1 2
        !   3 4
        ! and also a 2x1:
        !
        ! 1 3

        ASSERT(this%get_nsites() >= 4)

        ! set up the lattice indices, via the use of "k-vectors"
        temp_array(:, :) = 0
        if (this%t_bipartite_order) then
            if ( .not. (this%get_nsites() == 18 .or. this%get_nsites() == 8)) then
                call stop_all(this_routine, &
                    "bipartite only for 8 or 18 tilted sites for now")
            end if
            allocate(order(this%get_nsites()), source = 0)

            if (t_input_order) then
                order = orbital_order
            else
                if (this%get_nsites() == 18) then
                    order = [ 1, 2, 10, 3, 4, 11, 5, 12, 6, 13, 7, 14, 8, 15, 16, 9, 17, 18]
                else if (this%get_nsites() == 8) then
                    order = [1,2,5,3,6,4,7,8]
                end if
            end if
        else
            allocate(order(this%get_nsites()), source = [(i, i = 1, this%get_nsites())])
        end if


        k = 0
        l = 1
        m = this%get_nsites() / 2 + 1
        do i = -this%length(1) + 1, 0
            do j = -k, k

                temp_array(j, i) = order(l)
                l = l + 1

            end do
            k = k + 1
        end do

        ! here i need to change the k-vectors, differently, if it is a
        ! rectangular tilted lattice..
        ! and for now, until i have implemented it better over
        ! lattice vectors assume lx - ly <= 1! only one difference
        k = k - 1
        ! or should i do an inbetween-step if lx /= ly? this is also
        ! possible

        offset = abs(this%length(1) - this%length(2))
        k_min = -this%length(1) + 1
        k_max = this%length(2) - offset
        do i = 1, offset

            do j = k_min, k_max
                temp_array(j, i) = l
                l = l + 1
            end do

            ! shift the y indication by 1 up or down
            k_min = k_min + 1
            k_max = k_max + 1
        end do

        if (this%length(1) < this%length(2)) then
            k_min = k_min
            k_max = k_max - 1
        else if (this%length(1) > this%length(2)) then
            k_min = k_min + 1
            k_max = k_max
        else
            ! otherwise k_min and k_max where never defined
            k_min = -k
            k_max = k
        end if

        do i = offset + 1, this%length(2)


            ! if (this%t_bipartite_order) then
            !     do j = k_min, k_max, 2
            !
            !         temp_array(j, i) = m
            !
            !         m = m + 1
            !
            !     end do
            !     do j = k_min + 1, k_max - 1, 2
            !         temp_array(j, i) = l
            !         l = l + 1
            !     end do
            ! else
                do j = k_min, k_max

                    temp_array(j, i) = order(l)

                    l = l + 1

                end do
            ! end if

            ! k_min is negative
            k_min = k_min + 1
            k_max = k_max - 1
        end do

        up = cshift(temp_array, -1, 1)
        down = cshift(temp_array, 1, 1)
        right = cshift(temp_array, 1, 2)
        left = cshift(temp_array, -1, 2)

        ! apply the periodic boundary conditions to the neighbors
        pbc = this%length(1)
        ! for rectangular tilted lattices this is different of course

        if (this%length(1) == 1) then
            pbc_1 = [2, 0]
            pbc_2 = [this%length(2), -this%length(2)]
        else if (this%length(2) == 1) then
            pbc_1 = [this%length(1), this%length(1)]
            pbc_2 = [2, 0]
        else
            pbc_1 = [this%length(1), this%length(1)]
            pbc_2 = [this%length(2), -this%length(2)]
        end if

        ! do something like and do this generally maybe..
        call apply_pbc_tilted(up, pbc_1, pbc_2, up_ur, up_dr, up_ul, up_dl, up_rr, up_ll)
        call apply_pbc_tilted(down, pbc_1, pbc_2, down_ur, down_dr, down_ul, &
            down_dl, down_rr, down_ll)
        call apply_pbc_tilted(right, pbc_1, pbc_2, right_ur, right_dr, right_ul, &
            right_dl, right_rr, right_ll)
        call apply_pbc_tilted(left, pbc_1, pbc_2, left_ur, left_dr, left_ul, &
            left_dl, left_rr, left_ll)

        k = 0
        l = 1
        ! now get the neighbors
        if (this%is_periodic()) then
            ! fully periodic case
            do i = -this%length(1) + 1, 0
                do j = -k, k
                    ! make the neigbors list
                    up_nn = maxval([up(j, i), up_ur(j, i), up_dr(j, i), up_ul(j, i), &
                                    up_dl(j, i)])

                    if (up_nn == 0) then
                        up_nn = maxval([up_rr(j, i), up_ll(j, i)])
                        if (up_nn == 0) then
                            print *, " up: smth wrong!"
                        end if
                    end if

                    down_nn = maxval([down(j, i), down_ur(j, i), down_dr(j, i), &
                        down_ul(j, i), down_dl(j, i)])

                    if (down_nn == 0) then
                        down_nn = maxval([down_rr(j, i), down_ll(j, i)])

                        if (down_nn == 0) then
                            print *, "down: smth wrong!"
                        end if
                    end if

                    right_nn = maxval([right(j, i), right_ur(j, i), right_dr(j, i), &
                        right_ul(j, i), right_dl(j, i)])

                    if (right_nn == 0) then
                        right_nn = right_ll(j, i)
                        if (right_nn == 0) then
                            print *, "right: smth wrong!"
                        end if
                    end if

                    left_nn = maxval([left(j, i), left_ur(j, i), left_dr(j, i), &
                        left_ul(j, i), left_dl(j, i)])

                    if (left_nn == 0) then
                        left_nn = left_rr(j, i)
                        if (left_nn == 0) then
                            print *, "left: smth wrong!"
                        end if
                    end if

                    neigh = sort_unique([up_nn, down_nn, left_nn, right_nn])
                    ! also start to store the k-vector here!
                    ! have to be sure that i make it correct
                    k_vec = [i, j, 0]
                    r_vec = [j, i, 0]
                    this%sites(order(l)) = site(order(l), size(neigh), neigh, k_vec, r_vec)

                    l = l + 1

                    deallocate(neigh)

                end do
                k = k + 1
            end do

            k = k - 1

            k_min = -this%length(1) + 1
            k_max = this%length(2) - offset
            do i = 1, offset
                do j = k_min, k_max
                    ! make the neigbors list
                    up_nn = maxval([up(j, i), up_ur(j, i), up_dr(j, i), up_ul(j, i), &
                                    up_dl(j, i)])

                    if (up_nn == 0) then
                        up_nn = maxval([up_rr(j, i), up_ll(j, i)])
                        if (up_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    down_nn = maxval([down(j, i), down_ur(j, i), down_dr(j, i), &
                        down_ul(j, i), down_dl(j, i)])

                    if (down_nn == 0) then
                        down_nn = maxval([down_rr(j, i), down_ll(j, i)])

                        if (down_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    right_nn = maxval([right(j, i), right_ur(j, i), right_dr(j, i),&
                        right_ul(j, i), right_dl(j, i)])

                    if (right_nn == 0) then
                        right_nn = right_ll(j, i)
                        if (right_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    left_nn = maxval([left(j, i), left_ur(j, i), left_dr(j, i), &
                        left_ul(j, i), left_dl(j, i)])

                    if (left_nn == 0) then
                        left_nn = left_rr(j, i)
                        if (left_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    neigh = sort_unique([up_nn, down_nn, left_nn, right_nn])

                    k_vec = [i, j, 0]
                    r_vec = [j, 1, 0]
                    this%sites(order(l)) = site(order(l), size(neigh), neigh, k_vec, r_vec)

                    l = l + 1

                    deallocate(neigh)

                end do
                k_min = k_min + 1
                k_max = k_max + 1
            end do

            if (this%length(1) < this%length(2)) then
                k_min = k_min
                k_max = k_max - 1
            else if (this%length(1) > this%length(2)) then
                k_min = k_min + 1
                k_max = k_max
            else
                k_min = -k
                k_max = k
            end if

            do i = offset + 1, this%length(2)
                do j = k_min, k_max
                    ! make the neigbors list
                    up_nn = maxval([up(j, i), up_ur(j, i), up_dr(j, i), up_ul(j, i), &
                                    up_dl(j, i)])

                    if (up_nn == 0) then
                        up_nn = maxval([up_rr(j, i), up_ll(j, i)])
                        if (up_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    down_nn = maxval([down(j, i), down_ur(j, i), down_dr(j, i), &
                        down_ul(j, i), down_dl(j, i)])

                    if (down_nn == 0) then
                        down_nn = maxval([down_rr(j, i), down_ll(j, i)])

                        if (down_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    right_nn = maxval([right(j, i), right_ur(j, i), right_dr(j, i), &
                        right_ul(j, i), right_dl(j, i)])

                    if (right_nn == 0) then
                        right_nn = right_ll(j, i)
                        if (right_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    left_nn = maxval([left(j, i), left_ur(j, i), left_dr(j, i), &
                        left_ul(j, i), left_dl(j, i)])

                    if (left_nn == 0) then
                        left_nn = left_rr(j, i)
                        if (left_nn == 0) then
                            print *, "smth wrong!"
                        end if
                    end if

                    neigh = sort_unique([up_nn, down_nn, left_nn, right_nn])

                    k_vec = [i, j, 0]
                    r_vec = [j, i, 0]
                    this%sites(order(l)) = site(order(l), size(neigh), neigh, k_vec, r_vec)

                    l = l + 1

                    deallocate(neigh)

                end do
                k_min = k_min + 1
                k_max = k_max - 1

            end do
        else if (this%is_periodic(1)) then
            ! only apply (x,x) periodicity
            do i = -this%length(1) + 1, 0
                do j = -k, k

                    up_nn = maxval([up(j, i), up_ur(j, i), up_dl(j, i)])
                    down_nn = maxval([down(j, i), down_ur(j, i), down_dl(j, i)])
                    left_nn = maxval([left(j, i), left_ur(j, i), left_dl(j, i)])
                    right_nn = maxval([right(j, i), right_ur(j, i), right_dl(j, i)])

                    temp_neigh = [up_nn, down_nn, left_nn, right_nn]

                    neigh = sort_unique(pack(temp_neigh, temp_neigh > 0))

                    this%sites(order(l)) = site(order(l), size(neigh), neigh)

                    l = l + 1

                    deallocate(neigh)
                end do
                k = k + 1
            end do
            k = k - 1

            do i = 1, this%length(1)
                do j = -k, k

                    up_nn = maxval([up(j, i), up_ur(j, i), up_dl(j, i)])
                    down_nn = maxval([down(j, i), down_ur(j, i), down_dl(j, i)])
                    left_nn = maxval([left(j, i), left_ur(j, i), left_dl(j, i)])
                    right_nn = maxval([right(j, i), right_ur(j, i), right_dl(j, i)])

                    temp_neigh = [up_nn, down_nn, left_nn, right_nn]

                    neigh = sort_unique(pack(temp_neigh, temp_neigh > 0))

                    this%sites(order(l)) = site(order(l), size(neigh), neigh)

                    l = l + 1

                    deallocate(neigh)
                end do
                k = k - 1
            end do

        else if (this%is_periodic(2)) then
            ! only apply (x,-x) periodicity
            do i = -this%length(1) + 1, 0
                do j = -k, k

                    up_nn = maxval([up(j, i), up_ul(j, i), up_dr(j, i)])
                    down_nn = maxval([down(j, i), down_ul(j, i), down_dr(j, i)])
                    left_nn = maxval([left(j, i), left_ul(j, i), left_dr(j, i)])
                    right_nn = maxval([right(j, i), right_ul(j, i), right_dr(j, i)])

                    temp_neigh = [up_nn, down_nn, left_nn, right_nn]

                    neigh = sort_unique(pack(temp_neigh, temp_neigh > 0))

                    this%sites(order(l)) = site(order(l), size(neigh), neigh)

                    l = l + 1

                    deallocate(neigh)
                end do
                k = k + 1
            end do
            k = k - 1

            do i = 1, this%length(1)
                do j = -k, k

                    up_nn = maxval([up(j, i), up_ul(j, i), up_dr(j, i)])
                    down_nn = maxval([down(j, i), down_ul(j, i), down_dr(j, i)])
                    left_nn = maxval([left(j, i), left_ul(j, i), left_dr(j, i)])
                    right_nn = maxval([right(j, i), right_ul(j, i), right_dr(j, i)])

                    temp_neigh = [up_nn, down_nn, left_nn, right_nn]

                    neigh = sort_unique(pack(temp_neigh, temp_neigh > 0))

                    this%sites(order(l)) = site(order(l), size(neigh), neigh)

                    l = l + 1

                    deallocate(neigh)
                end do
                k = k - 1
            end do

        else
            ! non-periodic case
            do i = -this%length(1) + 1, 0
                do j = -k, k
                    ! only a neighbor if the index is non-zero!
                    temp_neigh = [up(j, i), down(j, i), left(j, i), right(j, i)]

                    neigh = sort_unique(pack(temp_neigh, temp_neigh > 0))

                    this%sites(order(l)) = site(order(l), size(neigh), neigh)

                    l = l + 1

                    deallocate(neigh)

                end do
                k = k + 1
            end do

            k = k - 1
            do i = 1, this%length(1)
                do j = -k, k
                    ! only a neighbor if the index is non-zero!
                    temp_neigh = [up(j, i), down(j, i), left(j, i), right(j, i)]

                    neigh = sort_unique(pack(temp_neigh, temp_neigh > 0))

                    this%sites(order(l)) = site(order(l), size(neigh), neigh)

                    l = l + 1

                    deallocate(neigh)

                end do
                k = k - 1
            end do
        end if

    end subroutine init_sites_tilted