1!--------------------------------------------------------------------------------------------------! 2! DFTB+: general package for performing fast atomistic simulations ! 3! Copyright (C) 2006 - 2019 DFTB+ developers group ! 4! ! 5! See the LICENSE file for terms of usage and distribution. ! 6!--------------------------------------------------------------------------------------------------! 7 8!> Contains a type with basis information 9module dftbp_orbitals 10 use dftbp_accuracy 11 use dftbp_message 12 use dftbp_constants, only : shellNames 13 implicit none 14 private 15 16 public :: TOrbitals, getShellNames, orbitalNames 17 18 19 !> Contains information about the orbitals of the species/atoms in the system 20 type TOrbitals 21 22 !> Nr. of shells for each atomic species (nSpecies) 23 integer, allocatable :: nShell(:) 24 25 !> Nr. of orbitals for each atomic species (nSpecies) 26 integer, allocatable :: nOrbSpecies(:) 27 28 !> Nr. of orbitals for each atom (nAtom) 29 integer, allocatable :: nOrbAtom(:) 30 31 !> Ang. momentum of the a particular l-shell on a particular species (maxval(nShell), nSpecies) 32 integer, allocatable :: angShell(:,:) 33 34 !> The shell which contains the given orbital on an atom 35 !> (maxval(nOrbSpecies), nSpecies) 36 integer, allocatable :: iShellOrb(:,:) 37 38 !> Starting pos. within the atomic block of the each of the shells of each species 39 !> (maxval(nShell)+1, nSpecies) 40 integer, allocatable :: posShell(:,:) 41 42 !> Max. nr. of shells for any species 43 integer :: mShell 44 45 !> Max. nr. of orbitals for any species 46 integer :: mOrb 47 48 !> Total number of orbitals in system. 49 integer :: nOrb 50 end type TOrbitals 51 52 !> Length of labels for atomic orbitals 53 integer, parameter :: lenOrbitalNames = 9 54 55 !> Names of the atomic orbitals in tesseral basis 56 !> general set for f orbitals (not cubic), see 57 !> http://winter.group.shef.ac.uk/orbitron/AOs/4f/equations.html 58 character(lenOrbitalNames), parameter :: orbitalNames(-3:3,0:3) = reshape([& 59 & ' ',' ',' ',' ',' ',' ',' ',& 60 & ' ',' ','y ','z ','x ',' ',' ',& 61 & ' ','xy ','yz ','z2 ','xz ','x2-y2 ',' ',& 62 & 'y(3x2-y2)','x2+y2+z2 ','yz2 ','z3 ','xz2 ','z(x2-y2) ','x(x2-3y2)'& 63 &], [7,4]) 64 65contains 66 67 !> Builds a unique names for the atomic orbitals 68 !> Assign 's', 'p', 'd' to first occurring shells, then 's2', 'p2', ... 69 subroutine getShellNames(iSpecie, orb, shellNamesTmp) 70 71 !> atomic specie 72 integer, intent(in) :: iSpecie 73 74 !> orbital info 75 type(TOrbitals), intent(in) :: orb 76 77 !> output string naming the atomic orbital 78 character(sc), intent(out), allocatable :: shellNamesTmp(:) 79 80 integer :: ii, jj 81 integer, allocatable :: ind(:) 82 character(sc) :: sindx 83 84 !allocate(names(orb%nShell(iSpecie))) 85 allocate(shellNamesTmp(orb%nShell(iSpecie))) 86 allocate(ind(orb%nShell(iSpecie))) 87 ind = 1 88 89 do ii = 1, orb%nShell(iSpecie) 90 write(shellNamesTmp(ii), "(A)") shellNames(orb%angShell(ii, iSpecie) + 1) 91 if (any(shellNamesTmp(ii) == shellNamesTmp(1:ii-1))) then 92 ! at least one example of this shell already 93 ind(ii) = ind(ii) + 1 94 write(sindx,'(I0)') ind(ii) 95 if (len(trim(adjustl(shellNamesTmp(ii)))) + len(trim(sindx)) > sc) then 96 call error("Shell labels are too long: "//trim(adjustl(shellNamesTmp(ii)))//trim(sindx)) 97 else 98 shellNamesTmp(ii) = trim(adjustl(shellNamesTmp(ii)))//trim(sindx) 99 end if 100 end if 101 end do 102 deallocate(ind) 103 104 end subroutine getShellNames 105 106 107end module dftbp_orbitals 108