TABLE OF CONTENTS
- ABINIT/m_wvl_descr_psp
- defs_wvltypes/wvl_descr_atoms_set
- defs_wvltypes/wvl_descr_atoms_set_sym
- defs_wvltypes/wvl_descr_free
- defs_wvltypes/wvl_descr_psp_fill
- defs_wvltypes/wvl_descr_psp_set
ABINIT/m_wvl_descr_psp [ Modules ]
NAME
wvl_descr_psp
FUNCTION
COPYRIGHT
Copyright (C) 2008-2024 ABINIT group (DC) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_wvl_descr_psp 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 28 use defs_datatypes, only : pseudopotential_type, pseudopotential_gth_type 29 30 implicit none 31 32 private
defs_wvltypes/wvl_descr_atoms_set [ Functions ]
[ Top ] [ defs_wvltypes ] [ Functions ]
NAME
wvl_descr_atoms_set
FUNCTION
Defines wvl%atoms% data structure
INPUTS
acell(3)=unit cell length scales (bohr) dtset <type(dataset_type)>=all input variables for this dataset
OUTPUT
wvl <type(wvl_internal_type)>= wavelet type | nat = number of atoms | ntypes = number of species | alat1 = acell(1) | alat2 = acell(2) | alat3 = acell(3) | iatype = types for atoms | lfrztyp = flag for the movement of atoms. | natpol = integer related to polarisation at the first step
SOURCE
361 subroutine wvl_descr_atoms_set(acell, icoulomb, natom, ntypat, typat, wvl) 362 363 use defs_wvltypes 364 #if defined HAVE_BIGDFT 365 use BigDFT_API, only: atoms_data_null,f_routine,f_release_routine,& 366 & astruct_set_n_atoms,astruct_set_n_types,& 367 & allocate_atoms_nat,allocate_atoms_ntypes 368 #endif 369 implicit none 370 371 !Arguments ------------------------------------ 372 !scalars 373 integer, intent(in) :: icoulomb, natom, ntypat 374 type(wvl_internal_type), intent(inout) :: wvl 375 !arrays 376 integer, intent(in) :: typat(natom) 377 real(dp), intent(in) :: acell(3) 378 379 !Local variables------------------------------- 380 !scalars 381 #if defined HAVE_BIGDFT 382 integer :: itype 383 #endif 384 385 ! ********************************************************************* 386 387 #if defined HAVE_BIGDFT 388 389 call f_routine('wvl_descr_atoms_set') 390 391 wvl%atoms=atoms_data_null() 392 393 !We create the atoms_data structure from this dataset 394 !to be used later in BigDFT routines. 395 if (icoulomb == 0) then 396 wvl%atoms%astruct%geocode = 'P' 397 else if (icoulomb == 1) then 398 wvl%atoms%astruct%geocode = 'F' 399 else if (icoulomb == 2) then 400 wvl%atoms%astruct%geocode = 'S' 401 end if 402 write(wvl%atoms%astruct%units, "(A)") "Bohr" 403 404 call astruct_set_n_atoms(wvl%atoms%astruct, natom) 405 call astruct_set_n_types(wvl%atoms%astruct, ntypat) 406 407 do itype = 1, ntypat, 1 408 write(wvl%atoms%astruct%atomnames(itype), "(A,I2)") "At. type", itype 409 end do 410 wvl%atoms%astruct%cell_dim(1) = acell(1) 411 wvl%atoms%astruct%cell_dim(2) = acell(2) 412 wvl%atoms%astruct%cell_dim(3) = acell(3) 413 wvl%atoms%astruct%iatype = typat 414 415 wvl%atoms%astruct%sym%symObj = 0 416 417 call allocate_atoms_nat(wvl%atoms) 418 call allocate_atoms_ntypes(wvl%atoms) 419 420 call f_release_routine() 421 422 #else 423 BIGDFT_NOTENABLED_ERROR() 424 if (.false.) write(std_out,*) icoulomb,natom,ntypat,wvl%h(1),typat(1),acell(1) 425 #endif 426 427 end subroutine wvl_descr_atoms_set
defs_wvltypes/wvl_descr_atoms_set_sym [ Functions ]
[ Top ] [ defs_wvltypes ] [ Functions ]
NAME
wvl_descr_atoms_set_sym
FUNCTION
Add symmetry information to wvl%atoms data structure.
INPUTS
acell(3)=unit cell length scales (bohr) dtset <type(dataset_type)>=all input variables for this dataset
OUTPUT
wvl <type(wvl_internal_type)>= wavelet type | nat = number of atoms | ntypes = number of species | alat1 = acell(1) | alat2 = acell(2) | alat3 = acell(3) | iatype = types for atoms | lfrztyp = flag for the movement of atoms. | natpol = integer related to polarisation at the first step
SOURCE
455 subroutine wvl_descr_atoms_set_sym(wvl, efield, irrzon, nsppol, nsym, phnons, & 456 & symafm, symrel, tnons, tolsym) 457 458 use defs_wvltypes 459 use m_ab7_symmetry 460 #if defined HAVE_BIGDFT 461 use BigDFT_API, only: astruct_set_symmetries 462 #endif 463 implicit none 464 465 !Arguments ------------------------------------ 466 !scalars 467 integer, intent(in) :: nsppol,nsym 468 real(dp), intent(in) :: tolsym 469 type(wvl_internal_type), intent(inout) :: wvl 470 !arrays 471 integer, intent(in) :: symafm(3,3,nsym), symrel(3,3,nsym) 472 integer, target, intent(in) :: irrzon(:,:,:) 473 real(dp), target, intent(in) :: phnons(:,:,:) 474 real(dp), intent(in) :: efield(3), tnons(3,nsym) 475 476 !Local variables------------------------------- 477 !scalars 478 #if defined HAVE_BIGDFT 479 integer :: errno 480 character(len=500) :: message 481 #endif 482 483 ! ********************************************************************* 484 485 #if defined HAVE_BIGDFT 486 487 write(message, '(a,a)' ) ch10,& 488 & ' wvl_descr_atoms_set_sym: Create symmetries for the wvl%atoms object.' 489 call wrtout(std_out,message,'COLL') 490 491 wvl%atoms%astruct%sym%symObj = -1 492 nullify(wvl%atoms%astruct%sym%irrzon) 493 nullify(wvl%atoms%astruct%sym%phnons) 494 call astruct_set_symmetries(wvl%atoms%astruct, (nsym <= 1), tolsym, efield, nsppol) 495 if (nsym > 1) then 496 call symmetry_set_n_sym(wvl%atoms%astruct%sym%symObj, nsym, symrel, tnons, symafm, errno) 497 end if 498 wvl%atoms%astruct%sym%irrzon => irrzon 499 wvl%atoms%astruct%sym%phnons => phnons 500 501 #else 502 BIGDFT_NOTENABLED_ERROR() 503 if (.false.) write(std_out,*) nsppol,nsym,tolsym,wvl%h(1),symafm(1,1,1),symrel(1,1,1),& 504 & irrzon(1,1,1),phnons(1,1,1),efield(1),tnons(1,1) 505 #endif 506 507 end subroutine wvl_descr_atoms_set_sym
defs_wvltypes/wvl_descr_free [ Functions ]
[ Top ] [ defs_wvltypes ] [ Functions ]
NAME
wvl_descr_free
FUNCTION
Free the wvl%atoms% datastructure (deallocate or nullify)
INPUTS
wvl <type(wvl_internal_type)>=internal variables for wavelets
OUTPUT
wvl <type(wvl_internal_type)>=internal variables for wavelets
SOURCE
302 subroutine wvl_descr_free(wvl) 303 304 use defs_wvltypes 305 #if defined HAVE_BIGDFT 306 use BigDFT_API, only : deallocate_atoms_data 307 use dynamic_memory 308 #endif 309 implicit none 310 311 !Arguments ------------------------------------ 312 !scalars 313 type(wvl_internal_type), intent(inout) :: wvl 314 !arrays 315 316 !Local variables------------------------------- 317 !scalars 318 319 ! ********************************************************************* 320 321 #if defined HAVE_BIGDFT 322 !These arrays are pointers on memory handled by ABINIT. 323 nullify(wvl%atoms%astruct%sym%irrzon) 324 nullify(wvl%atoms%astruct%sym%phnons) 325 if (associated(wvl%atoms%nlccpar)) then 326 call f_free_ptr(wvl%atoms%nlccpar) 327 end if 328 call deallocate_atoms_data(wvl%atoms) 329 #endif 330 if(allocated(wvl%npspcode_paw_init_guess)) then 331 ABI_FREE(wvl%npspcode_paw_init_guess) 332 end if 333 end subroutine wvl_descr_free
defs_wvltypes/wvl_descr_psp_fill [ Functions ]
[ Top ] [ defs_wvltypes ] [ Functions ]
NAME
wvl_descr_psp_fill
FUNCTION
Read the radii from @pspunit, if any and fill values with default otherwise.
INPUTS
OUTPUT
SOURCE
161 subroutine wvl_descr_psp_fill(gth_params, ipsp, ixc, nelpsp, nzatom, pspunit) 162 163 #if defined HAVE_BIGDFT 164 use BigDFT_API, only: atomic_info, UNINITIALIZED, psp_from_data 165 #endif 166 implicit none 167 168 !Arguments ------------------------------------ 169 integer, intent(in) :: ipsp, pspunit, nzatom, nelpsp, ixc 170 type(pseudopotential_gth_type), intent(inout) :: gth_params 171 !Local variables------------------------------- 172 #if defined HAVE_BIGDFT 173 integer :: ios, ii, nzatom_, nelpsp_, npspcode_, ixc_ 174 real(dp) :: ehomo, radfine 175 logical :: exists 176 character(len = 2) :: symbol 177 character(len=100) :: line 178 character(len=500) :: message 179 #endif 180 181 ! *************************************************************************** 182 183 #if defined HAVE_BIGDFT 184 185 ! check if gth_params%psppar have been set 186 if (any(gth_params%psppar == UNINITIALIZED(1._dp))) then 187 call atomic_info(nzatom, nelpsp, symbol = symbol) 188 ixc_ = ixc 189 if (ixc_>0.and.ixc_< 10) ixc_=1 190 if (ixc_>0.and.ixc_>=10) ixc_=11 191 call psp_from_data(symbol, nzatom_, nelpsp_, npspcode_, ixc_, & 192 & gth_params%psppar(:,:,ipsp), exists) 193 if(.not. exists) then 194 write(message,'(a,a,a,a)')ch10,& 195 & "wvl_descr_psp_fill : bug, chemical element not found in BigDFT table",ch10,& 196 & "Action: upgrade BigDFT table" 197 call wrtout(ab_out,message,'COLL') 198 ABI_BUG(message) 199 end if 200 gth_params%set(ipsp) = .true. 201 end if 202 203 ! Try to read radii from pspunit 204 if (pspunit /= 0) then 205 read (pspunit, '(a100)', iostat = ios) line 206 if (ios /= 0) then 207 line='' 208 end if 209 ! We try to read the values from the pseudo. 210 read (line, *, iostat = ios) gth_params%radii_cf(ipsp, 1), gth_params%radii_cf(ipsp, 2), & 211 & gth_params%radii_cf(ipsp, 3) 212 if (ios /= 0 .or. gth_params%radii_cf(ipsp, 3) < zero) then 213 read (line, *, iostat = ios) gth_params%radii_cf(ipsp, 1), gth_params%radii_cf(ipsp, 2) 214 gth_params%radii_cf(ipsp, 3) = UNINITIALIZED(gth_params%radii_cf(ipsp, 3)) 215 end if 216 if (ios /= 0 .or. gth_params%radii_cf(ipsp, 1) < zero .or. gth_params%radii_cf(ipsp, 2) < zero) then 217 gth_params%radii_cf(ipsp, 1) = UNINITIALIZED(gth_params%radii_cf(ipsp, 1)) 218 gth_params%radii_cf(ipsp, 2) = UNINITIALIZED(gth_params%radii_cf(ipsp, 2)) 219 write(message, '(a,a,a,a,a,a,a)' ) '-', ch10,& 220 & '- wvl_descr_psp_fill : COMMENT -',ch10,& 221 & "- the pseudo-potential does not include geometric information,",ch10,& 222 & '- values have been computed.' 223 call wrtout(ab_out,message,'COLL') 224 call wrtout(std_out, message,'COLL') 225 end if 226 end if 227 228 ! Update radii. 229 if (gth_params%radii_cf(ipsp, 1) == UNINITIALIZED(gth_params%radii_cf(ipsp, 1))) then 230 call atomic_info(nzatom, nelpsp, ehomo = ehomo) 231 !assigning the radii by calculating physical parameters 232 gth_params%radii_cf(ipsp, 1)=1._dp/sqrt(abs(2._dp*ehomo)) 233 end if 234 if (gth_params%radii_cf(ipsp, 2) == UNINITIALIZED(gth_params%radii_cf(ipsp, 2))) then 235 radfine = gth_params%psppar(0, 0, ipsp) 236 do ii = 1, 4, 1 237 if (gth_params%psppar(ii, 0, ipsp) /= zero) then 238 radfine = min(radfine, gth_params%psppar(ii, 0, ipsp)) 239 end if 240 end do 241 gth_params%radii_cf(ipsp, 2)=radfine 242 end if 243 if (gth_params%radii_cf(ipsp, 3) == UNINITIALIZED(gth_params%radii_cf(ipsp, 3))) then 244 gth_params%radii_cf(ipsp, 3)=gth_params%radii_cf(ipsp, 2) 245 end if 246 247 ! Set flag. 248 if (gth_params%radii_cf(ipsp, 1) >= 0.d0 .and. gth_params%radii_cf(ipsp, 2) >= 0.d0) then 249 gth_params%hasGeometry(ipsp) = .true. 250 else 251 gth_params%hasGeometry(ipsp) = .false. 252 end if 253 254 write(message, '(a,f12.7,a,f12.7,a,f12.7)' )& 255 & ' radii_cf(1)=', gth_params%radii_cf(ipsp, 1),& 256 & '; radii_cf(2)=', gth_params%radii_cf(ipsp, 2),& 257 & '; rad_cov=', gth_params%radii_cf(ipsp, 3) 258 call wrtout(ab_out,message,'COLL') 259 call wrtout(std_out, message,'COLL') 260 261 ! Some consistency checks on radii, to be moved earlier, but need to update test refs. 262 ! gth_params%radii_cf(ipsp, 3) = min(gth_params%radii_cf(ipsp, 3), gth_params%radii_cf(ipsp, 2)) 263 264 ! This was one before 265 ! maxrad=zero 266 ! do ii=0,2,1 267 ! if (ii==1) maxrad=zero 268 ! if (gth_params%psppar(ii,0,ipsp)/=zero) maxrad=max(maxrad,gth_params%psppar(ii,0,ipsp)) 269 ! end do 270 ! if (maxrad== zero) then 271 ! gth_params%radii_cf(ipsp,3)=zero 272 ! else 273 ! gth_params%radii_cf(ipsp,3)=max( & 274 !& min(dtset%wvl_crmult*psps%gth_params%radii_cf(ipsp,1), & 275 !& 15._dp*maxrad)/dtset%wvl_frmult,gth_params%radii_cf(ipsp,2)) 276 ! end if 277 278 #else 279 BIGDFT_NOTENABLED_ERROR() 280 if (.false.) write(std_out,*) ipsp,pspunit,nzatom,nelpsp,ixc,gth_params%psppar 281 #endif 282 283 end subroutine wvl_descr_psp_fill
defs_wvltypes/wvl_descr_psp_set [ Functions ]
[ Top ] [ defs_wvltypes ] [ Functions ]
NAME
wvl_descr_psp_set
FUNCTION
Defines the part of the wvl%atoms%-datastructure which depends on psps
INPUTS
psps <type(pseudopotential_type)>=variables related to pseudopotentials dtset <type(dataset_type)>=all input variables for this dataset
OUTPUT
nsppol=number of spin components spinat(3,natom)=spin polarization on each atom wvl <type(wvl_internal_type)> = wavelet type | psppar = The covalence radii for each pseudo | pspcod = the format -or code- of psp generation | iasctype = semicore code (see defs_datatype) | nzatom = charge of the nucleus | nelpsp = the ionic pseudo-charge | natsc = number of atoms with semicore
SOURCE
69 subroutine wvl_descr_psp_set(filoccup, nsppol, psps, spinat, wvl) 70 71 use defs_wvltypes 72 #if defined HAVE_BIGDFT 73 use BigDFT_API, only: aoig_set,UNINITIALIZED,dict_init,dict_free,dictionary, & 74 & operator(//), bigdft_mpi, dict_set 75 use BigDFT_API, only: psp_data_merge_to_dict, psp_dict_fill_all, atomic_info, & 76 & psp_dict_analyse, atomic_data_set_from_dict, merge_input_file_to_dict 77 #endif 78 implicit none 79 80 !Arguments ------------------------------------ 81 !scalars 82 integer, intent(in) :: nsppol 83 type(wvl_internal_type), intent(inout) :: wvl 84 type(pseudopotential_type), intent(in) :: psps 85 character(len = *), intent(in) :: filoccup 86 !arrays 87 real(dp),intent(in) :: spinat(:,:) 88 89 !Local variables------------------------------- 90 #if defined HAVE_BIGDFT 91 integer :: ityp,pspcod 92 logical :: exists 93 real(dp) :: radii_cf(3) 94 character(len=2) :: symbol 95 character(len=27) :: filename 96 type(dictionary), pointer :: dict 97 #endif 98 99 ! ********************************************************************* 100 101 #if defined HAVE_BIGDFT 102 103 !We create the atoms_data structure, the part that is dependent from psp. 104 do ityp=1,size(psps%pspcod) 105 wvl%atoms%psppar(:,:,ityp)= psps%gth_params%psppar(:,:,ityp) 106 wvl%atoms%npspcode(ityp) = merge(psps%pspcod(ityp),7,psps%pspcod(ityp)/=17) 107 wvl%atoms%ixcpsp(ityp) = psps%pspxc(ityp) 108 wvl%atoms%nzatom(ityp) = int(psps%znucltypat(ityp)) 109 wvl%atoms%nelpsp(ityp) = int(psps%ziontypat(ityp)) 110 end do 111 wvl%atoms%donlcc = .false. 112 113 call dict_init(dict) 114 radii_cf = UNINITIALIZED(1._dp) 115 do ityp = 1, wvl%atoms%astruct%ntypes, 1 116 call atomic_info(wvl%atoms%nzatom(ityp), wvl%atoms%nelpsp(ityp), & 117 & symbol = symbol) 118 write(wvl%atoms%astruct%atomnames(ityp), "(A)") symbol 119 filename = 'psppar.' // trim(symbol) 120 pspcod=psps%pspcod(ityp);if (pspcod==17) pspcod=7 ! PAW specificity 121 call psp_data_merge_to_dict(dict // filename, int(psps%znucltypat(ityp)), & 122 & int(psps%ziontypat(ityp)), pspcod, psps%pspxc(ityp), & 123 & psps%gth_params%psppar(0:4,0:6,ityp), radii_cf, & 124 & UNINITIALIZED(1._dp), UNINITIALIZED(1._dp)) 125 call psp_dict_fill_all(dict, trim(symbol), psps%pspxc(ityp)) 126 end do 127 128 call psp_dict_analyse(dict, wvl%atoms) 129 inquire(file = filoccup, exist = exists) 130 if (exists) then 131 call merge_input_file_to_dict(dict//"Atomic occupation",filoccup,bigdft_mpi) 132 end if 133 ! Need to add spinat in the dictionary, later 134 if (spinat(1,1) == zero .and. .false.) then 135 call dict_set(dict // "spinat", spinat(1,1)) 136 end if 137 call atomic_data_set_from_dict(dict, "Atomic occupation", wvl%atoms, nsppol) 138 call dict_free(dict) 139 140 #else 141 BIGDFT_NOTENABLED_ERROR() 142 if (.false.) write(std_out,*) nsppol,wvl%h(1),psps%npsp,filoccup,spinat(1,1) 143 #endif 144 145 end subroutine wvl_descr_psp_set