1! 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt. 6! See Docs/Contributors.txt for a list of contributors. 7! 8 module chemical 9 10 use sys 11 use precision, only: dp 12 13 implicit none 14 15 private 16 17 public :: atomic_number 18 public :: number_of_species, species_label 19 public :: is_floating, is_bessel, is_synthetic 20 public :: read_chemical_types, print_chemical_type 21 22 ! Public due to Bcast routines 23 public :: chemical_types, chemical_list 24 25 ! Species information 26 type chemical_types 27 integer :: no_of_species 28 character(len=20), pointer :: spec_label(:) 29 integer, pointer :: z(:) 30 end type chemical_types 31 32 type(chemical_types), save :: chemical_list 33 34 35 CONTAINS 36 37 subroutine check(i) 38 integer, intent(in) :: i 39 if (i.lt.0 .or. i.gt.chemical_list%no_of_species) 40 $ call die("Wrong species number requested") 41 end subroutine check 42 43 function number_of_species() 44 integer number_of_species 45 number_of_species = chemical_list%no_of_species 46 end function number_of_species 47 48 function species_label(i) 49 character(len=20) species_label 50 integer, intent(in) :: i 51 52 call check(i) 53 species_label = chemical_list%spec_label(i) 54 end function species_label 55 56 function atomic_number(i) 57 integer atomic_number 58 integer, intent(in) :: i 59 60 call check(i) 61 atomic_number = chemical_list%z(i) 62 end function atomic_number 63! ------- 64 function is_floating(i) 65 logical is_floating 66 integer, intent(in) :: i 67 68 call check(i) 69 is_floating = (chemical_list%z(i) .le. 0) 70 end function is_floating 71! ------- 72 73 function is_bessel(i) 74 logical is_bessel 75 integer, intent(in) :: i 76 77 call check(i) 78 is_bessel = (chemical_list%z(i) .eq. -100) 79 end function is_bessel 80! ------- 81! Checks whether we are dealing with a synthetic atom 82! 83 function is_synthetic(i) 84 logical is_synthetic 85 integer, intent(in) :: i 86 87 call check(i) 88 ! Note that we could have a synthetic ghost atom, with 89 ! z <= -200 90 is_synthetic = (abs(chemical_list%z(i)) .gt. 200) 91 end function is_synthetic 92!--- 93 subroutine read_chemical_types(silent) 94 95 use parallel, only : Node 96 use fdf 97 98 logical, intent(in), optional :: silent ! default .false. 99 100 integer nsp, isp 101 integer ns_read 102 103 type(block_fdf) :: bfdf 104 type(parsed_line), pointer :: pline 105 106 character(len=20) :: label 107 character(len=256) :: msg 108 109 integer :: z, is 110 logical :: found 111 logical :: lsilent 112 113 ! Determine whether we should be silent 114 lsilent = .false. 115 if ( present(silent) ) lsilent = silent 116 if ( Node /= 0 ) lsilent = .true. 117 118 ! Default to 0 119 nsp = fdf_integer('Number_of_species',0) 120 121 ! The most important thing to find is 122 ! the block containing the species 123 found = fdf_block('Chemical_species_label', bfdf) 124 if (.not. found ) 125 $ call die("Block Chemical_species_label does not exist.") 126 127 if ( nsp == 0 ) then 128 ns_read = fdf_block_linecount('Chemical_species_label', 'iin') 129 else 130 ns_read = nsp 131 end if 132 ! If they are not equal we notify the user 133 if ( nsp /= ns_read ) then 134 nsp = ns_read 135 end if 136 if ( nsp == 0 ) call die("No species found!!!") 137 138 allocate(chemical_list%spec_label(nsp)) 139 allocate(chemical_list%z(nsp)) 140 chemical_list%no_of_species = nsp 141 142 ns_read = 0 143 do while( fdf_bline(bfdf,pline) ) 144 if ( .not. fdf_bmatch(pline,'iin') ) cycle 145 146 ns_read = ns_read + 1 147 148 ! Get species information 149 isp = fdf_bintegers(pline,1) 150 label = fdf_bnames(pline,1) 151 z = fdf_bintegers(pline,2) 152 153 ! We cannot test label names in this 154 ! loop as isp may be non-linear 155 if ( isp < 1 .or. nsp < isp ) 156 $ call die("Wrong specnum in Chemical_species_label") 157 158 chemical_list%z(isp) = z 159 chemical_list%spec_label(isp) = label 160 161 end do 162 if ( ns_read /= nsp ) 163 & call die("Not enough species in block") 164 165 if ( .not. lsilent ) then 166 ! Align output, always 167 do isp = 1 , nsp 168 call print_chemical_type(isp) 169 end do 170 write(*,*) ! new-line 171 end if 172 173 ! Check that none of the chemical species are the 174 ! same 175 if ( nsp > 1 ) then 176 do z = 1 , nsp - 1 177 do is = z+1, nsp 178 if (trim(species_label(z))==trim(species_label(is))) then 179 write(msg,'(2(a,i0,2a),2a)') 180 & "Specie index/label = ", 181 & z,'/',trim(species_label(z)), " has same label as ", 182 & is,'/',trim(species_label(is)),". ", 183 & " Use a different one for hygienic reasons." 184 call die(trim(msg)) 185 end if 186 end do 187 end do 188 end if 189 190 end subroutine read_chemical_types 191 192 subroutine print_chemical_type(isp) 193 integer, intent(in) :: isp 194 195 character(len=256) :: label 196 integer :: z 197 logical :: floating, bessel, synthetic 198 199 label = species_label(isp) 200 z = atomic_number(isp) 201 202 floating = z <= 0 203 bessel = z == -100 204 synthetic = abs(z) > 200 205 206 ! Bessel *must* be checked first! 207 if ( bessel ) then 208 write(6,'(a,i3,3a)') 'Species number: ',isp, 209 $ ' Label: ', trim(label), 210 $ ' (floating Bessel functions)' 211 else if ( floating ) then 212 write(6,'(a,i3,a,i4,3a)') 'Species number: ',isp, 213 $ ' Atomic number: ', z, 214 $ ' Label: ', trim(label), 215 $ ' (floating PAOs)' 216 else 217 write(6,'(a,i3,a,i4,2a)') 'Species number: ',isp, 218 $ ' Atomic number: ', z, 219 $ ' Label: ', trim(label) 220 end if 221 end subroutine print_chemical_type 222 223 end module chemical 224