TABLE OF CONTENTS
ABINIT/m_wvl_projectors [ Modules ]
NAME
m_wvl_projectors
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_projectors 23 24 use defs_basis 25 use defs_wvltypes 26 use m_errors 27 use m_abicore 28 use m_atomdata 29 30 use defs_datatypes, only : pseudopotential_type 31 use m_geometry, only : xred2xcart 32 33 implicit none 34 35 private
ABINIT/wvl_projectors_free [ Functions ]
NAME
wvl_projectors_free
FUNCTION
Freeing routine.
INPUTS
OUTPUT
SIDE EFFECTS
proj <type(wvl_projectors_type)>=projectors information in a wavelet basis.
SOURCE
169 subroutine wvl_projectors_free(proj) 170 171 #if defined HAVE_BIGDFT 172 use BigDFT_API, only : free_DFT_PSP_projectors,deallocate_gwf 173 #endif 174 implicit none 175 176 !Arguments ------------------------------------ 177 !scalars 178 type(wvl_projectors_type),intent(inout) :: proj 179 180 !Local variables ------------------------- 181 #if defined HAVE_BIGDFT 182 integer :: ii 183 #endif 184 185 ! ********************************************************************* 186 187 #if defined HAVE_BIGDFT 188 189 call free_DFT_PSP_projectors(proj%nlpsp) 190 191 if (allocated(proj%G)) then 192 do ii=1,size(proj%G) 193 ! MT dec 2014: cannot call bigdft deallocation routine 194 ! because content of proj%G datastructure was created 195 ! without f_malloc (without memory profiling). 196 ! call deallocate_gwf(proj%G(ii)) 197 if (associated(proj%G(ii)%ndoc)) then 198 ABI_FREE(proj%G(ii)%ndoc) 199 end if 200 if (associated(proj%G(ii)%nam)) then 201 ABI_FREE(proj%G(ii)%nam) 202 end if 203 if (associated(proj%G(ii)%nshell)) then 204 ABI_FREE(proj%G(ii)%nshell) 205 end if 206 if (associated(proj%G(ii)%psiat)) then 207 ABI_FREE(proj%G(ii)%psiat) 208 end if 209 if (associated(proj%G(ii)%xp)) then 210 ABI_FREE(proj%G(ii)%xp) 211 end if 212 end do 213 ABI_FREE(proj%G) 214 end if 215 216 #else 217 BIGDFT_NOTENABLED_ERROR() 218 if (.false.) write(std_out,*) proj%nlpsp 219 #endif 220 221 end subroutine wvl_projectors_free
ABINIT/wvl_projectors_set [ Functions ]
NAME
wvl_projectors_set
FUNCTION
Allocate and compute the access keys for the projectors when the positions of the atoms are given. The array to store projectors is also allocated, use wvl_projectors_free() to free them after use.
INPUTS
dtset <type(dataset_type)>=internal variables used by wavelets, describing | wvl_internal=desciption of the wavelet box. | natom=number of atoms. mpi_enreg=information about MPI parallelization psps <type(pseudopotential_type)>=variables related to pseudopotentials rprimd(3,3)=dimensional primitive translations in real space (bohr) xred(3,natom)=reduced dimensionless atomic coordinates
OUTPUT
proj <type(wvl_projector_type)>=projectors information for wavelets. | keys=its access keys for compact storage.
SIDE EFFECTS
SOURCE
72 subroutine wvl_projectors_set(me, natom, proj, psps, rprimd, wfs, wvl, wvl_frmult, xred) 73 74 #if defined HAVE_BIGDFT 75 use BigDFT_API, only: createProjectorsArrays, wvl_timing => timing 76 #endif 77 implicit none 78 79 !Arguments ------------------------------------ 80 !scalars 81 integer, intent(in) :: natom, me 82 real(dp), intent(in) :: wvl_frmult 83 type(pseudopotential_type),intent(in) :: psps 84 type(wvl_projectors_type),intent(inout) :: proj 85 type(wvl_wf_type),intent(in) :: wfs 86 type(wvl_internal_type), intent(in) :: wvl 87 !arrays 88 real(dp),intent(in) :: rprimd(3,3),xred(3,natom) 89 90 !Local variables------------------------------- 91 !scalars 92 #if defined HAVE_BIGDFT 93 integer :: idata 94 logical,parameter :: wvl_debug=.false. 95 character(len=500) :: message 96 !arrays 97 real(dp),allocatable :: xcart(:,:) 98 #endif 99 100 ! ********************************************************************* 101 102 #if defined HAVE_BIGDFT 103 !Consistency checks, are all pseudo true GTH pseudo with geometric information? 104 do idata = 1, psps%npsp, 1 105 if (.not. psps%gth_params%set(idata)) then 106 write(message, '(a,a,a,a,I0,a,a,a)' ) ch10,& 107 & ' wvl_projectors_set : consistency checks failed,', ch10, & 108 & ' no GTH parameters found for type number ', idata, '.', ch10, & 109 & ' Check your input pseudo files.' 110 ABI_ERROR(message) 111 end if 112 if (.not. psps%gth_params%hasGeometry(idata)) then 113 write(message, '(a,a,a,a,a,a)' ) ch10,& 114 & ' wvl_projectors_set : consistency checks failed,', ch10, & 115 & ' the given GTH parameters has no geometry information.', ch10, & 116 & ' Upgrade your input pseudo files to GTH with geometric informatoins.' 117 ABI_ERROR(message) 118 end if 119 end do 120 121 if (wvl_debug) then 122 call wvl_timing(me,'CrtProjectors ','ON') 123 end if 124 125 !Store xcart for each atom 126 ABI_MALLOC(xcart,(3, natom)) 127 call xred2xcart(natom, rprimd, xcart, xred) 128 call createProjectorsArrays(wfs%ks%Lzd%Glr,xcart,wvl%atoms,wfs%ks%orbs,& 129 psps%gth_params%radii_cf,wvl_frmult,wvl_frmult,wvl%h(1),wvl%h(2),& 130 wvl%h(3),.false.,proj%nlpsp,proj%G) 131 write(message, '(a,a,a,a,I0)' ) ch10,& 132 & ' wvl_projectors_set : allocate projectors data,', ch10, & 133 & ' size of the compressed array: ', proj%nlpsp%nprojel 134 call wrtout(std_out,message,'COLL') 135 136 !Deallocations 137 ABI_FREE(xcart) 138 139 if (wvl_debug) then 140 call wvl_timing(me,'CrtProjectors ','OF') 141 end if 142 143 #else 144 BIGDFT_NOTENABLED_ERROR() 145 if (.false.) write(std_out,*) natom,me,wvl_frmult,psps%npsp,proj%nlpsp,wfs%ks,wvl%h(1),& 146 & rprimd(1,1),xred(1,1) 147 #endif 148 149 end subroutine wvl_projectors_set