TABLE OF CONTENTS
ABINIT/ddb_interpolate [ Functions ]
NAME
ddb_interpolate
FUNCTION
Interpolate the ddb onto a fine set of q-points using the interatomic force constants and write the result in a DDB file.
INPUTS
OUTPUT
NOTES
SOURCE
71 subroutine ddb_interpolate(ifc, crystal, inp, ddb, ddb_hdr, asrq0, prefix, comm) 72 73 !Arguments ------------------------------- 74 !scalars 75 type(ifc_type),intent(in) :: ifc 76 type(anaddb_dataset_type),target,intent(inout) :: inp 77 type(crystal_t),intent(in) :: crystal 78 type(ddb_type),intent(inout) :: ddb 79 type(ddb_hdr_type),intent(inout) :: ddb_hdr 80 type(asrq0_t),intent(inout) :: asrq0 81 integer,intent(in) :: comm 82 character(len=*),intent(in) :: prefix 83 !arrays 84 85 !Local variables ------------------------- 86 !scalars 87 integer,parameter :: master=0 88 integer :: nsym,natom,ntypat,mband,nqpt_fine 89 integer :: msize,nsize,mpert,nblok,mtyp 90 integer :: rftyp 91 integer :: ii,iblok,jblok,iqpt,ipert1,ipert2,idir1,idir2 92 integer :: nprocs,my_rank 93 character(len=500) :: msg 94 character(len=fnlen) :: ddb_out_filename, ddb_out_nc_filename 95 type(ddb_type) :: ddb_new 96 !arrays 97 integer :: rfphon(4),rfelfd(4),rfstrs(4) 98 integer,allocatable :: blkflg(:,:,:,:) 99 real(dp) :: qpt(3), qptnrm(3), qpt_padded(3,3) 100 real(dp),allocatable :: d2cart(:,:,:,:,:),d2red(:,:,:,:,:) 101 real(dp),pointer :: qpt_fine(:,:) 102 103 ! ********************************************************************* 104 105 106 ! Only master works for the time being 107 nprocs = xmpi_comm_size(comm); my_rank = xmpi_comm_rank(comm) 108 if (my_rank /= master) return 109 110 ! =================================== 111 ! Copy dimensions and allocate arrays 112 ! =================================== 113 114 write(msg, '(a,(80a),a,a,a,a)' ) ch10,('=',ii=1,80),ch10,ch10,' Treat the first list of vectors ',ch10 115 call wrtout(std_out,msg,'COLL') 116 call wrtout(ab_out,msg,'COLL') 117 118 nullify(qpt_fine) 119 nqpt_fine = inp%nph1l 120 qpt_fine => inp%qph1l 121 122 rftyp=inp%rfmeth 123 124 nsym = Crystal%nsym 125 natom = Crystal%natom 126 ntypat = Crystal%ntypat 127 128 mband = ddb_hdr%mband 129 130 mtyp = max(ddb_hdr%mblktyp, 2) ! Limited to 2nd derivatives of total energy 131 ddb_hdr%mblktyp = mtyp 132 133 mpert = ddb%mpert 134 msize = 3 * mpert * 3 * mpert !; if (mtyp==3) msize=msize*3*mpert 135 nsize = 3 * mpert * 3 * mpert 136 nblok = nqpt_fine 137 138 ddb_new%nblok = nblok 139 call ddb_new%malloc(msize,nblok,natom,ntypat,mpert) 140 ddb_new%flg = 0 141 ddb_new%amu = ddb%amu 142 ddb_new%typ = 1 143 ddb_new%qpt = zero 144 ddb_new%nrm = one 145 146 ABI_MALLOC(d2cart,(2,3,mpert,3,mpert)) 147 ABI_MALLOC(d2red,(2,3,mpert,3,mpert)) 148 ABI_MALLOC(blkflg,(3,mpert,3,mpert)) 149 150 blkflg = 1 151 152 rfphon(1:2)=1; rfelfd(1:2)=0; rfstrs(1:2)=0 153 qpt_padded = zero 154 155 156 ddb_hdr%dscrpt = 'Interpolated DDB using interatomic force constants' 157 ddb_hdr%nblok = nblok 158 159 ! ================================================ 160 ! Interpolate the dynamical matrix at each q-point 161 ! ================================================ 162 163 qptnrm = one 164 165 do iqpt=1,nqpt_fine 166 167 ! Initialisation of the phonon wavevector 168 qpt(:)=qpt_fine(:,iqpt) 169 170 if (inp%nph1l /= 0) qptnrm(1) = inp%qnrml1(iqpt) 171 172 ! Look for the information in the DDB 173 qpt_padded(:,1) = qpt 174 call ddb%get_block(iblok,qpt_padded,qptnrm,rfphon,rfelfd,rfstrs,rftyp) 175 176 if (iblok /= 0) then 177 178 ! q-point is present in the ddb. No interpolation needed. 179 180 d2cart(1,:,:,:,:) = reshape(ddb%val(1,:,iblok), shape = (/3,mpert,3,mpert/)) 181 d2cart(2,:,:,:,:) = reshape(ddb%val(2,:,iblok), shape = (/3,mpert,3,mpert/)) 182 else 183 184 ! Get d2cart using the interatomic forces and the 185 ! long-range coulomb interaction through Ewald summation 186 call gtdyn9(ddb%acell,Ifc%atmfrc,Ifc%dielt,Ifc%dipdip,Ifc%dyewq0,d2cart, & 187 crystal%gmet,ddb%gprim,ddb%mpert,natom,Ifc%nrpt,qptnrm(1), & 188 qpt, crystal%rmet,ddb%rprim,Ifc%rpt,Ifc%trans,crystal%ucvol, & 189 Ifc%wghatm,crystal%xred,ifc%zeff,ifc%qdrp_cart,ifc%ewald_option,xmpi_comm_self, & 190 dipquad=Ifc%dipquad,quadquad=Ifc%quadquad) 191 192 end if 193 194 ! Eventually impose the acoustic sum rule based on previously calculated d2asr 195 call asrq0%apply(natom, ddb%mpert, ddb%msize, crystal%xcart, d2cart) 196 197 ! Transform d2cart into reduced coordinates. 198 call d2cart_to_red(d2cart,d2red,crystal%gprimd,crystal%rprimd,mpert, & 199 & natom,ntypat,crystal%typat,crystal%ucvol,crystal%zion) 200 201 ! Store the dynamical matrix into a block of the new ddb 202 jblok = iqpt 203 ddb_new%val(1,1:nsize,jblok) = reshape(d2red(1,:,:,:,:), shape = (/nsize/)) 204 ddb_new%val(2,1:nsize,jblok) = reshape(d2red(2,:,:,:,:), shape = (/nsize/)) 205 206 ! Store the q-point 207 ddb_new%qpt(1:3,jblok) = qpt 208 ddb_new%nrm(1,jblok) = qptnrm(1) 209 210 ! Set the flags 211 ii=0 212 do ipert2=1,mpert 213 do idir2=1,3 214 do ipert1=1,mpert 215 do idir1=1,3 216 ii=ii+1 217 if (ipert1<=natom.and.ipert2<=natom) then 218 ddb_new%flg(ii,jblok) = 1 219 end if 220 end do 221 end do 222 end do 223 end do 224 225 end do ! iqpt 226 227 ! Copy the flags for Gamma 228 qpt_padded(:,1) = zero 229 qptnrm = one 230 call ddb%get_block(iblok,qpt_padded,qptnrm,rfphon,rfelfd,rfstrs,rftyp) 231 call ddb_new%get_block(jblok,qpt_padded,qptnrm,rfphon,rfelfd,rfstrs,rftyp) 232 233 if (iblok /= 0 .and. jblok /= 0) then 234 235 ii=0 236 do ipert2=1,mpert 237 do idir2=1,3 238 do ipert1=1,mpert 239 do idir1=1,3 240 ii=ii+1 241 ddb_new%flg(ii,jblok) = ddb%flg(ii,iblok) 242 end do 243 end do 244 end do 245 end do 246 247 end if 248 249 250 ! ======================= 251 ! Write the new DDB files 252 ! ======================= 253 254 if (my_rank == master) then 255 256 ddb_out_filename = strcat(prefix, "_DDB") 257 258 call ddb_new%write_txt(ddb_hdr, ddb_out_filename) 259 260 ! Write one separate nc file for each q-point 261 do jblok=1,nblok 262 write(ddb_out_nc_filename,'(2a,i5.5,a)') trim(prefix),'_qpt_',jblok,'_DDB.nc' 263 call ddb_new%write_nc(ddb_hdr, ddb_out_nc_filename, jblok) 264 end do 265 266 end if 267 268 ! =========== 269 ! Free memory 270 ! =========== 271 272 call ddb_new%free() 273 ABI_FREE(d2cart) 274 ABI_FREE(d2red) 275 ABI_FREE(blkflg) 276 277 end subroutine ddb_interpolate
ABINIT/m_ddb_interpolate [ Modules ]
NAME
m_ddb_interpolate
FUNCTION
Interpolate the ddb onto a fine set of q-points using the interatomic force constants and write the result in a DDB file.
COPYRIGHT
Copyright (C) 2008-2022 ABINIT group (GA) 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
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 module m_ddb_interpolate 24 25 use defs_basis 26 use m_errors 27 use m_xmpi 28 use m_abicore 29 use m_ddb 30 use m_ddb_hdr 31 use m_ifc 32 use m_nctk 33 #ifdef HAVE_NETCDF 34 use netcdf 35 #endif 36 37 use m_anaddb_dataset, only : anaddb_dataset_type 38 use m_crystal, only : crystal_t 39 use m_io_tools, only : get_unit 40 use m_fstrings, only : strcat 41 use m_dynmat, only : gtdyn9, d2cart_to_red 42 43 implicit none 44 45 private