TABLE OF CONTENTS
- ABINIT/m_primcell_ddb_info
- m_primcell_ddb_info/destroy_primcell_ddb_info
- m_primcell_ddb_info/init_primcell_ddb_info
- m_primcell_ddb_info/read_primcell_ddb_info
- m_primcell_ddb_info/write_primcell_ddb_info
ABINIT/m_primcell_ddb_info [ Modules ]
NAME
m_primcell_ddb_info
FUNCTION
Module for container object passing cell information to SC phonon calculation in abinit Container type is defined, and destruction
COPYRIGHT
Copyright (C) 1999-2024 ABINIT group (MJV) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
SOURCE
20 #if defined HAVE_CONFIG_H 21 #include "config.h" 22 #endif 23 24 #include "abi_common.h" 25 26 27 module m_primcell_ddb_info 28 29 use defs_basis 30 use m_abicore 31 use m_errors 32 33 use m_io_tools, only : open_file 34 35 implicit none 36 37 type primcell_ddb_info 38 ! scalars 39 integer :: brav,mpert,msym,natom,nrpt,nsym,ntypat 40 integer :: dipdip 41 real(dp) :: ucvol 42 43 ! arrays 44 integer , allocatable :: indsym(:,:,:) ! indsym(4,nsym,natom)=label given by subroutine symatm 45 integer , allocatable :: symrec(:,:,:) ! (3,3,nsym) recip space symops 46 integer , allocatable :: symrel(:,:,:) ! (3,3,nsym) real space symops 47 integer , allocatable :: typat(:) ! typat(natom)=integer label of each type of atom (1,2,...) 48 49 real(dp), allocatable :: acell(:) ! acell(3) 50 real(dp), allocatable :: amu(:) ! amu(ntypat)=mass of the atoms (atomic mass unit) 51 real(dp), allocatable :: dielt(:,:) ! dielt(3,3)=dielectric tensor 52 real(dp), allocatable :: dyewq0(:,:,:) ! dyewq0(3,3,natom)=Ewald part of the dynamical matrix, at q=0 53 real(dp), allocatable :: gmet(:,:) ! gmet(3,3) 54 real(dp), allocatable :: gprim(:,:) ! gprim(3,3) 55 real(dp), allocatable :: rcan(:,:) ! rcan(3,natom)=atomic position in canonical coordinates 56 real(dp), allocatable :: rmet(:,:) ! rmet(3,3) 57 real(dp), allocatable :: rprim(:,:) ! rprim(3,3)=dimensionless primitive translations in real space 58 real(dp), allocatable :: rpt(:,:) ! rpt(3,nrpt)=canonical coordinates of the R points in the unit cell 59 real(dp), allocatable :: trans(:,:) ! trans(3,natom)=atomic translations : xred = rcan + trans 60 real(dp), allocatable :: wghatm(:,:,:) ! wghatm(natom,natom,nrpt)=weights assoc to a pair of atoms + R vector 61 real(dp), allocatable :: xred(:,:) ! xred(3,natom)= relative coords of atoms in unit cell (dimensionless) 62 real(dp), allocatable :: zeff(:,:,:) ! zeff(3,3,natom)=effective charge on each atom 63 64 end type primcell_ddb_info 65 66 ! Now the subroutines in the module 67 contains
m_primcell_ddb_info/destroy_primcell_ddb_info [ Functions ]
[ Top ] [ m_primcell_ddb_info ] [ Functions ]
NAME
destroy_primcell_ddb_info
FUNCTION
deallocate stuoff in primcell_ddb_info
COPYRIGHT
Copyright (C) 1999-2024 ABINIT group (MJV) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
OUTPUT
NOTES
SOURCE
437 subroutine destroy_primcell_ddb_info (pcell) 438 439 use defs_basis 440 441 !Arguments ------------------------------------ 442 type(primcell_ddb_info), intent(inout) :: pcell 443 444 ! ************************************************************************* 445 if (allocated(pcell%indsym)) then 446 ABI_FREE(pcell%indsym) 447 end if 448 if (allocated(pcell%symrec)) then 449 ABI_FREE(pcell%symrec) 450 end if 451 if (allocated(pcell%symrel)) then 452 ABI_FREE(pcell%symrel) 453 end if 454 if (allocated(pcell%typat )) then 455 ABI_FREE(pcell%typat) 456 end if 457 if (allocated(pcell%acell )) then 458 ABI_FREE(pcell%acell) 459 end if 460 if (allocated(pcell%amu )) then 461 ABI_FREE(pcell%amu) 462 end if 463 if (allocated(pcell%dielt )) then 464 ABI_FREE(pcell%dielt) 465 end if 466 if (allocated(pcell%dyewq0)) then 467 ABI_FREE(pcell%dyewq0) 468 end if 469 if (allocated(pcell%gmet )) then 470 ABI_FREE(pcell%gmet) 471 end if 472 if (allocated(pcell%gprim )) then 473 ABI_FREE(pcell%gprim) 474 end if 475 if (allocated(pcell%rcan )) then 476 ABI_FREE(pcell%rcan) 477 end if 478 if (allocated(pcell%rmet )) then 479 ABI_FREE(pcell%rmet) 480 end if 481 if (allocated(pcell%rprim )) then 482 ABI_FREE(pcell%rprim) 483 end if 484 if (allocated(pcell%rpt )) then 485 ABI_FREE(pcell%rpt) 486 end if 487 if (allocated(pcell%trans )) then 488 ABI_FREE(pcell%trans) 489 end if 490 if (allocated(pcell%wghatm)) then 491 ABI_FREE(pcell%wghatm) 492 end if 493 if (allocated(pcell%xred )) then 494 ABI_FREE(pcell%xred) 495 end if 496 if (allocated(pcell%zeff )) then 497 ABI_FREE(pcell%zeff) 498 end if 499 500 end subroutine destroy_primcell_ddb_info 501 502 end module m_primcell_ddb_info
m_primcell_ddb_info/init_primcell_ddb_info [ Functions ]
[ Top ] [ m_primcell_ddb_info ] [ Functions ]
NAME
init_primcell_ddb_info
FUNCTION
init and fill primcell_ddb_info
COPYRIGHT
Copyright (C) 1999-2024 ABINIT group (MJV) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
OUTPUT
pcell= structure to allocate and fill
NOTES
SOURCE
93 subroutine init_primcell_ddb_info (pcell,brav,dipdip,mpert,msym,natom,nrpt,nsym,ntypat,ucvol,& 94 & indsym,symrec,symrel,typat,& 95 & acell,amu,dielt,dyewq0,gmet,gprim,rcan,rmet,rprim,rpt,trans,wghatm,xred,zeff) 96 97 use defs_basis 98 99 !Arguments ------------------------------------ 100 type(primcell_ddb_info), intent(inout) :: pcell 101 102 integer, intent(in) :: brav,mpert,msym,natom,nrpt,nsym,ntypat,dipdip 103 real(dp), intent(in) :: ucvol 104 105 integer, intent(in) :: indsym(4,nsym,natom) 106 integer, intent(in) :: symrec(3,3,nsym) 107 integer, intent(in) :: symrel(3,3,nsym) 108 integer, intent(in) :: typat(natom) 109 110 real(dp), intent(in) :: acell(3) 111 real(dp), intent(in) :: amu(ntypat) 112 real(dp), intent(in) :: dielt(3,3) 113 real(dp), intent(in) :: dyewq0(3,3,natom) 114 real(dp), intent(in) :: gmet(3,3) 115 real(dp), intent(in) :: gprim(3,3) 116 real(dp), intent(in) :: rcan(3,natom) 117 real(dp), intent(in) :: rmet(3,3) 118 real(dp), intent(in) :: rprim(3,3) 119 real(dp), intent(in) :: rpt(3,nrpt) 120 real(dp), intent(in) :: trans(3,natom) 121 real(dp), intent(in) :: wghatm(natom,natom,nrpt) 122 real(dp), intent(in) :: xred(3,natom) 123 real(dp), intent(in) :: zeff(3,3,natom) 124 125 !Local variables------------------------------- 126 127 ! ************************************************************************* 128 129 ! init dimensions 130 pcell%brav = brav 131 pcell%mpert = mpert 132 pcell%msym = msym 133 pcell%natom = natom 134 pcell%nrpt = nrpt 135 pcell%nsym = nsym 136 pcell%ntypat = ntypat 137 138 ! init scalar reals 139 pcell%ucvol = ucvol 140 141 ! allocate int 142 ABI_MALLOC(pcell%indsym,(4,nsym,natom)) 143 ABI_MALLOC(pcell%symrec,(3,3,nsym)) 144 ABI_MALLOC(pcell%symrel,(3,3,nsym)) 145 ABI_MALLOC(pcell%typat,(natom)) 146 147 ! allocate real 148 ABI_MALLOC(pcell%acell,(3)) 149 ABI_MALLOC(pcell%amu,(ntypat)) 150 ABI_MALLOC(pcell%dielt,(3,3)) 151 ABI_MALLOC(pcell%dyewq0,(3,3,natom)) 152 ABI_MALLOC(pcell%gmet,(3,3)) 153 ABI_MALLOC(pcell%gprim,(3,3)) 154 ABI_MALLOC(pcell%rcan,(3,natom)) 155 ABI_MALLOC(pcell%rmet,(3,3)) 156 ABI_MALLOC(pcell%rprim,(3,3)) 157 ABI_MALLOC(pcell%rpt,(3,nrpt)) 158 ABI_MALLOC(pcell%trans,(3,natom)) 159 ABI_MALLOC(pcell%wghatm,(natom,natom,nrpt)) 160 ABI_MALLOC(pcell%xred,(3,natom)) 161 ABI_MALLOC(pcell%zeff,(3,3,natom)) 162 163 ! init int 164 pcell%indsym = indsym 165 pcell%symrec = symrec 166 pcell%symrel = symrel 167 pcell%typat = typat 168 pcell%dipdip = dipdip 169 170 ! init real 171 pcell%acell = acell 172 pcell%amu = amu 173 pcell%dielt = dielt 174 pcell%dyewq0 = dyewq0 175 pcell%gmet = gmet 176 pcell%gprim = gprim 177 pcell%rcan = rcan 178 pcell%rmet = rmet 179 pcell%rprim = rprim 180 pcell%rpt = rpt 181 pcell%trans = trans 182 pcell%wghatm = wghatm 183 pcell%xred = xred 184 pcell%zeff = zeff 185 186 end subroutine init_primcell_ddb_info
m_primcell_ddb_info/read_primcell_ddb_info [ Functions ]
[ Top ] [ m_primcell_ddb_info ] [ Functions ]
NAME
read_primcell_ddb_info
FUNCTION
read in and fill primcell_ddb_info from the file name given in input
COPYRIGHT
Copyright (C) 1999-2024 ABINIT group (MJV) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
filename= name of file to read in
OUTPUT
t_primcell_ddb_info= structure to allocate and fill
NOTES
SOURCE
213 subroutine read_primcell_ddb_info (filename,pcell) 214 215 use defs_basis 216 217 !Arguments ------------------------------------ 218 character(len=*), intent(in) :: filename 219 type(primcell_ddb_info), intent(inout) :: pcell 220 221 !Local variables------------------------------- 222 integer :: unit 223 character(len=500) :: msg 224 character(len=13):: buffer 225 226 ! ************************************************************************* 227 228 if (open_file(filename,msg,newunit=unit) /= 0) then 229 ABI_ERROR(msg) 230 end if 231 232 ! read in dimensions 233 read(unit,'(a,I6)') buffer, pcell%brav 234 read(unit,'(a,I6)') buffer, pcell%dipdip 235 read(unit,'(a,I6)') buffer, pcell%mpert 236 read(unit,'(a,I6)') buffer, pcell%msym 237 read(unit,'(a,I6)') buffer, pcell%natom 238 read(unit,'(a,I6)') buffer, pcell%nrpt 239 read(unit,'(a,I6)') buffer, pcell%nsym 240 read(unit,'(a,I6)') buffer, pcell%ntypat 241 242 ! read in scalar reals 243 read(unit,'(a,E20.10)') buffer, pcell%ucvol 244 245 ! allocate int 246 ABI_MALLOC(pcell%indsym,(4,pcell%nsym,pcell%natom)) 247 ABI_MALLOC(pcell%symrec,(3,3,pcell%nsym)) 248 ABI_MALLOC(pcell%symrel,(3,3,pcell%nsym)) 249 ABI_MALLOC(pcell%typat,(pcell%natom)) 250 251 ! allocate real 252 ABI_MALLOC(pcell%acell,(3)) 253 ABI_MALLOC(pcell%amu,(pcell%ntypat)) 254 ABI_MALLOC(pcell%dielt,(3,3)) 255 ABI_MALLOC(pcell%dyewq0,(3,3,pcell%natom)) 256 ABI_MALLOC(pcell%gmet,(3,3)) 257 ABI_MALLOC(pcell%gprim,(3,3)) 258 ABI_MALLOC(pcell%rcan,(3,pcell%natom)) 259 ABI_MALLOC(pcell%rmet,(3,3)) 260 ABI_MALLOC(pcell%rprim,(3,3)) 261 ABI_MALLOC(pcell%rpt,(3,pcell%nrpt)) 262 ABI_MALLOC(pcell%trans,(3,pcell%natom)) 263 ABI_MALLOC(pcell%wghatm,(pcell%natom,pcell%natom,pcell%nrpt)) 264 ABI_MALLOC(pcell%xred,(3,pcell%natom)) 265 ABI_MALLOC(pcell%zeff,(3,3,pcell%natom)) 266 267 ! read in int 268 read(unit,'(a)') 269 read(unit,'(8I6)') pcell%indsym 270 read(unit,'(a)') 271 read(unit,'(8I6)') pcell%symrec 272 read(unit,'(a)') 273 read(unit,'(8I6)') pcell%symrel 274 read(unit,'(a)') 275 read(unit,'(8I6)') pcell%typat 276 277 ! read in real 278 read(unit,'(a)') 279 read(unit,'(3E20.10)') pcell%acell 280 read(unit,'(a)') 281 read(unit,'(3E20.10)') pcell%amu 282 read(unit,'(a)') 283 read(unit,'(3E20.10)') pcell%dielt 284 read(unit,'(a)') 285 read(unit,'(3E20.10)') pcell%dyewq0 286 read(unit,'(a)') 287 read(unit,'(3E20.10)') pcell%gmet 288 read(unit,'(a)') 289 read(unit,'(3E20.10)') pcell%gprim 290 read(unit,'(a)') 291 read(unit,'(3E20.10)') pcell%rcan 292 read(unit,'(a)') 293 read(unit,'(3E20.10)') pcell%rmet 294 read(unit,'(a)') 295 read(unit,'(3E20.10)') pcell%rprim 296 read(unit,'(a)') 297 read(unit,'(3E20.10)') pcell%rpt 298 read(unit,'(a)') 299 read(unit,'(3E20.10)') pcell%trans 300 read(unit,'(a)') 301 read(unit,'(3E20.10)') pcell%wghatm 302 read(unit,'(a)') 303 read(unit,'(3E20.10)') pcell%xred 304 read(unit,'(a)') 305 read(unit,'(3E20.10)') pcell%zeff 306 307 close(unit) 308 309 end subroutine read_primcell_ddb_info
m_primcell_ddb_info/write_primcell_ddb_info [ Functions ]
[ Top ] [ m_primcell_ddb_info ] [ Functions ]
NAME
write_primcell_ddb_info
FUNCTION
write out primcell_ddb_info to the file name given in input
COPYRIGHT
Copyright (C) 1999-2024 ABINIT group (MJV) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
filename= name of file to read in t_primcell_ddb_info= structure to allocate and fill
OUTPUT
writes to file
NOTES
SOURCE
338 subroutine write_primcell_ddb_info (filename,pcell) 339 340 use defs_basis 341 342 !Arguments ------------------------------------ 343 344 character(len=*), intent(in) :: filename 345 type(primcell_ddb_info), intent(in) :: pcell 346 347 !Local variables------------------------------- 348 integer :: unit 349 character(len=500) :: msg 350 351 ! ************************************************************************* 352 353 if (open_file(filename,msg,newunit=unit, form="formatted",status="unknown") /= 0) then 354 ABI_ERROR(msg) 355 end if 356 357 ! read out dimensions 358 write(unit,'(a,I6)') 'pcell%brav ', pcell%brav 359 write(unit,'(a,I6)') 'pcell%dipdip ', pcell%dipdip 360 write(unit,'(a,I6)') 'pcell%mpert ', pcell%mpert 361 write(unit,'(a,I6)') 'pcell%msym ', pcell%msym 362 write(unit,'(a,I6)') 'pcell%natom ', pcell%natom 363 write(unit,'(a,I6)') 'pcell%nrpt ', pcell%nrpt 364 write(unit,'(a,I6)') 'pcell%nsym ', pcell%nsym 365 write(unit,'(a,I6)') 'pcell%ntypat ', pcell%ntypat 366 367 ! write out scalar reals 368 write(unit,'(a,E20.10)') 'pcell%ucvol ', pcell%ucvol 369 370 ! write out int 371 write(unit,'(a)') 'pcell%indsym' 372 write(unit,'(8I6)') pcell%indsym 373 write(unit,'(a)') 'pcell%symrec' 374 write(unit,'(8I6)') pcell%symrec 375 write(unit,'(a)') 'pcell%symrel' 376 write(unit,'(8I6)') pcell%symrel 377 write(unit,'(a)') 'pcell%typat' 378 write(unit,'(8I6)') pcell%typat 379 380 ! write out real 381 write(unit,'(a)') 'pcell%acell' 382 write(unit,'(3E20.10)') pcell%acell 383 write(unit,'(a)') 'pcell%amu' 384 write(unit,'(3E20.10)') pcell%amu 385 write(unit,'(a)') 'pcell%dielt' 386 write(unit,'(3E20.10)') pcell%dielt 387 write(unit,'(a)') 'pcell%dyewq0' 388 write(unit,'(3E20.10)') pcell%dyewq0 389 write(unit,'(a)') 'pcell%gmet' 390 write(unit,'(3E20.10)') pcell%gmet 391 write(unit,'(a)') 'pcell%gprim' 392 write(unit,'(3E20.10)') pcell%gprim 393 write(unit,'(a)') 'pcell%rcan' 394 write(unit,'(3E20.10)') pcell%rcan 395 write(unit,'(a)') 'pcell%rmet' 396 write(unit,'(3E20.10)') pcell%rmet 397 write(unit,'(a)') 'pcell%rprim' 398 write(unit,'(3E20.10)') pcell%rprim 399 write(unit,'(a)') 'pcell%rpt' 400 write(unit,'(3E20.10)') pcell%rpt 401 write(unit,'(a)') 'pcell%trans' 402 write(unit,'(3E20.10)') pcell%trans 403 write(unit,'(a)') 'pcell%wghatm' 404 write(unit,'(3E20.10)') pcell%wghatm 405 write(unit,'(a)') 'pcell%xred' 406 write(unit,'(3E20.10)') pcell%xred 407 write(unit,'(a)') 'pcell%zeff' 408 write(unit,'(3E20.10)') pcell%zeff 409 410 close(unit) 411 412 end subroutine write_primcell_ddb_info