TABLE OF CONTENTS


ABINIT/m_primcell_ddb_info [ Modules ]

[ Top ] [ 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-2018 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

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 
28 module m_primcell_ddb_info
29 
30  use defs_basis
31  use m_abicore
32  use m_errors
33 
34  use m_io_tools, only : open_file
35 
36  implicit none
37 
38 type primcell_ddb_info
39 ! scalars
40   integer :: brav,mpert,msym,natom,nrpt,nsym,ntypat
41   integer :: dipdip
42   real(dp) :: ucvol
43 
44 ! arrays
45   integer , allocatable :: indsym(:,:,:) ! indsym(4,nsym,natom)=label given by subroutine symatm
46   integer , allocatable :: symrec(:,:,:) ! (3,3,nsym) recip space symops
47   integer , allocatable :: symrel(:,:,:) ! (3,3,nsym) real  space symops
48   integer , allocatable :: typat(:)      ! typat(natom)=integer label of each type of atom (1,2,...)
49 
50   real(dp), allocatable :: acell(:)      ! acell(3)
51   real(dp), allocatable :: amu(:)        ! amu(ntypat)=mass of the atoms (atomic mass unit)
52   real(dp), allocatable :: dielt(:,:)    ! dielt(3,3)=dielectric tensor
53   real(dp), allocatable :: dyewq0(:,:,:) ! dyewq0(3,3,natom)=Ewald part of the dynamical matrix, at q=0
54   real(dp), allocatable :: gmet(:,:)     ! gmet(3,3)
55   real(dp), allocatable :: gprim(:,:)    ! gprim(3,3)
56   real(dp), allocatable :: rcan(:,:)     ! rcan(3,natom)=atomic position in canonical coordinates
57   real(dp), allocatable :: rmet(:,:)     ! rmet(3,3)
58   real(dp), allocatable :: rprim(:,:)    ! rprim(3,3)=dimensionless primitive translations in real space
59   real(dp), allocatable :: rpt(:,:)      ! rpt(3,nrpt)=canonical coordinates of the R points in the unit cell
60   real(dp), allocatable :: trans(:,:)    ! trans(3,natom)=atomic translations : xred = rcan + trans
61   real(dp), allocatable :: wghatm(:,:,:) ! wghatm(natom,natom,nrpt)=weights assoc to a pair of atoms + R vector
62   real(dp), allocatable :: xred(:,:)     ! xred(3,natom)= relative coords of atoms in unit cell (dimensionless)
63   real(dp), allocatable :: zeff(:,:,:)   ! zeff(3,3,natom)=effective charge on each atom
64 
65 end type primcell_ddb_info
66 
67 ! Now the subroutines in the module
68 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-2018 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

PARENTS

CHILDREN

SOURCE

478 subroutine destroy_primcell_ddb_info (pcell)
479  
480  use defs_basis
481 
482 !This section has been created automatically by the script Abilint (TD).
483 !Do not modify the following lines by hand.
484 #undef ABI_FUNC
485 #define ABI_FUNC 'destroy_primcell_ddb_info'
486 !End of the abilint section
487 
488  implicit none
489 
490 !Arguments ------------------------------------
491  type(primcell_ddb_info), intent(inout) :: pcell
492 
493 ! *************************************************************************
494   if (allocated(pcell%indsym))  then
495     ABI_DEALLOCATE(pcell%indsym)
496   end if
497   if (allocated(pcell%symrec))  then
498     ABI_DEALLOCATE(pcell%symrec)
499   end if
500   if (allocated(pcell%symrel))  then
501     ABI_DEALLOCATE(pcell%symrel)
502   end if
503   if (allocated(pcell%typat ))  then
504     ABI_DEALLOCATE(pcell%typat)
505   end if
506   if (allocated(pcell%acell ))  then
507     ABI_DEALLOCATE(pcell%acell)
508   end if
509   if (allocated(pcell%amu   ))  then
510     ABI_DEALLOCATE(pcell%amu)
511   end if
512   if (allocated(pcell%dielt ))  then
513     ABI_DEALLOCATE(pcell%dielt)
514   end if
515   if (allocated(pcell%dyewq0))  then
516     ABI_DEALLOCATE(pcell%dyewq0)
517   end if
518   if (allocated(pcell%gmet  ))  then
519     ABI_DEALLOCATE(pcell%gmet)
520   end if
521   if (allocated(pcell%gprim ))  then
522     ABI_DEALLOCATE(pcell%gprim)
523   end if
524   if (allocated(pcell%rcan  ))  then
525     ABI_DEALLOCATE(pcell%rcan)
526   end if
527   if (allocated(pcell%rmet  ))  then
528     ABI_DEALLOCATE(pcell%rmet)
529   end if
530   if (allocated(pcell%rprim ))  then
531     ABI_DEALLOCATE(pcell%rprim)
532   end if
533   if (allocated(pcell%rpt   ))  then
534     ABI_DEALLOCATE(pcell%rpt)
535   end if
536   if (allocated(pcell%trans ))  then
537     ABI_DEALLOCATE(pcell%trans)
538   end if
539   if (allocated(pcell%wghatm))  then
540     ABI_DEALLOCATE(pcell%wghatm)
541   end if
542   if (allocated(pcell%xred  ))  then
543     ABI_DEALLOCATE(pcell%xred)
544   end if
545   if (allocated(pcell%zeff  ))  then
546     ABI_DEALLOCATE(pcell%zeff)
547   end if
548 
549  end subroutine destroy_primcell_ddb_info
550 
551 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-2018 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

PARENTS

CHILDREN

SOURCE

 98 subroutine init_primcell_ddb_info (pcell,brav,dipdip,mpert,msym,natom,nrpt,nsym,ntypat,ucvol,&
 99 &    indsym,symrec,symrel,typat,&
100 &    acell,amu,dielt,dyewq0,gmet,gprim,rcan,rmet,rprim,rpt,trans,wghatm,xred,zeff)
101  
102  use defs_basis
103 
104 !This section has been created automatically by the script Abilint (TD).
105 !Do not modify the following lines by hand.
106 #undef ABI_FUNC
107 #define ABI_FUNC 'init_primcell_ddb_info'
108 !End of the abilint section
109 
110  implicit none
111 
112 !Arguments ------------------------------------
113  type(primcell_ddb_info), intent(inout) :: pcell
114 
115  integer, intent(in) :: brav,mpert,msym,natom,nrpt,nsym,ntypat,dipdip
116  real(dp), intent(in) :: ucvol
117 
118  integer, intent(in) :: indsym(4,nsym,natom)
119  integer, intent(in) :: symrec(3,3,nsym)
120  integer, intent(in) :: symrel(3,3,nsym)
121  integer, intent(in) :: typat(natom)
122 
123  real(dp), intent(in) :: acell(3)
124  real(dp), intent(in) :: amu(ntypat)
125  real(dp), intent(in) :: dielt(3,3)
126  real(dp), intent(in) :: dyewq0(3,3,natom)
127  real(dp), intent(in) :: gmet(3,3)
128  real(dp), intent(in) :: gprim(3,3)
129  real(dp), intent(in) :: rcan(3,natom)
130  real(dp), intent(in) :: rmet(3,3)
131  real(dp), intent(in) :: rprim(3,3)
132  real(dp), intent(in) :: rpt(3,nrpt)
133  real(dp), intent(in) :: trans(3,natom)
134  real(dp), intent(in) :: wghatm(natom,natom,nrpt)
135  real(dp), intent(in) :: xred(3,natom)
136  real(dp), intent(in) :: zeff(3,3,natom)
137 
138 !Local variables-------------------------------
139 
140 ! *************************************************************************
141 
142 ! init dimensions
143   pcell%brav = brav
144   pcell%mpert = mpert
145   pcell%msym = msym
146   pcell%natom = natom
147   pcell%nrpt = nrpt
148   pcell%nsym = nsym
149   pcell%ntypat = ntypat
150 
151 ! init scalar reals
152   pcell%ucvol = ucvol
153 
154 ! allocate int
155   ABI_ALLOCATE(pcell%indsym,(4,nsym,natom))
156   ABI_ALLOCATE(pcell%symrec,(3,3,nsym))
157   ABI_ALLOCATE(pcell%symrel,(3,3,nsym))
158   ABI_ALLOCATE(pcell%typat,(natom))
159 
160 ! allocate real
161   ABI_ALLOCATE(pcell%acell,(3))
162   ABI_ALLOCATE(pcell%amu,(ntypat))
163   ABI_ALLOCATE(pcell%dielt,(3,3))
164   ABI_ALLOCATE(pcell%dyewq0,(3,3,natom))
165   ABI_ALLOCATE(pcell%gmet,(3,3))
166   ABI_ALLOCATE(pcell%gprim,(3,3))
167   ABI_ALLOCATE(pcell%rcan,(3,natom))
168   ABI_ALLOCATE(pcell%rmet,(3,3))
169   ABI_ALLOCATE(pcell%rprim,(3,3))
170   ABI_ALLOCATE(pcell%rpt,(3,nrpt))
171   ABI_ALLOCATE(pcell%trans,(3,natom))
172   ABI_ALLOCATE(pcell%wghatm,(natom,natom,nrpt))
173   ABI_ALLOCATE(pcell%xred,(3,natom))
174   ABI_ALLOCATE(pcell%zeff,(3,3,natom))
175 
176 ! init int
177   pcell%indsym = indsym
178   pcell%symrec = symrec
179   pcell%symrel = symrel
180   pcell%typat = typat
181   pcell%dipdip = dipdip
182 
183 ! init real
184   pcell%acell = acell
185   pcell%amu = amu
186   pcell%dielt = dielt
187   pcell%dyewq0 = dyewq0
188   pcell%gmet = gmet
189   pcell%gprim = gprim
190   pcell%rcan = rcan
191   pcell%rmet = rmet
192   pcell%rprim = rprim
193   pcell%rpt = rpt
194   pcell%trans = trans
195   pcell%wghatm = wghatm
196   pcell%xred = xred
197   pcell%zeff = zeff
198   
199  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-2018 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

PARENTS

CHILDREN

SOURCE

230 subroutine read_primcell_ddb_info (filename,pcell)
231  
232  use defs_basis
233 
234 !This section has been created automatically by the script Abilint (TD).
235 !Do not modify the following lines by hand.
236 #undef ABI_FUNC
237 #define ABI_FUNC 'read_primcell_ddb_info'
238 !End of the abilint section
239 
240  implicit none
241 
242 !Arguments ------------------------------------
243  character(len=*), intent(in) :: filename
244  type(primcell_ddb_info), intent(inout) :: pcell
245 
246 !Local variables-------------------------------
247  integer :: unit
248  character(len=500) :: msg
249  character(len=13):: buffer
250  
251 ! *************************************************************************
252 
253   if (open_file(filename,msg,newunit=unit) /= 0) then
254     MSG_ERROR(msg)
255   end if
256 
257 ! read in dimensions
258   read(unit,'(a,I6)') buffer, pcell%brav
259   read(unit,'(a,I6)') buffer, pcell%dipdip
260   read(unit,'(a,I6)') buffer, pcell%mpert
261   read(unit,'(a,I6)') buffer, pcell%msym
262   read(unit,'(a,I6)') buffer, pcell%natom
263   read(unit,'(a,I6)') buffer, pcell%nrpt
264   read(unit,'(a,I6)') buffer, pcell%nsym
265   read(unit,'(a,I6)') buffer, pcell%ntypat
266 
267 ! read in scalar reals
268   read(unit,'(a,E20.10)') buffer, pcell%ucvol
269 
270 ! allocate int
271   ABI_ALLOCATE(pcell%indsym,(4,pcell%nsym,pcell%natom))
272   ABI_ALLOCATE(pcell%symrec,(3,3,pcell%nsym))
273   ABI_ALLOCATE(pcell%symrel,(3,3,pcell%nsym))
274   ABI_ALLOCATE(pcell%typat,(pcell%natom))
275 
276 ! allocate real
277   ABI_ALLOCATE(pcell%acell,(3))
278   ABI_ALLOCATE(pcell%amu,(pcell%ntypat))
279   ABI_ALLOCATE(pcell%dielt,(3,3))
280   ABI_ALLOCATE(pcell%dyewq0,(3,3,pcell%natom))
281   ABI_ALLOCATE(pcell%gmet,(3,3))
282   ABI_ALLOCATE(pcell%gprim,(3,3))
283   ABI_ALLOCATE(pcell%rcan,(3,pcell%natom))
284   ABI_ALLOCATE(pcell%rmet,(3,3))
285   ABI_ALLOCATE(pcell%rprim,(3,3))
286   ABI_ALLOCATE(pcell%rpt,(3,pcell%nrpt))
287   ABI_ALLOCATE(pcell%trans,(3,pcell%natom))
288   ABI_ALLOCATE(pcell%wghatm,(pcell%natom,pcell%natom,pcell%nrpt))
289   ABI_ALLOCATE(pcell%xred,(3,pcell%natom))
290   ABI_ALLOCATE(pcell%zeff,(3,3,pcell%natom))
291 
292 ! read in int
293   read(unit,'(a)')
294   read(unit,'(8I6)') pcell%indsym
295   read(unit,'(a)')
296   read(unit,'(8I6)') pcell%symrec
297   read(unit,'(a)')
298   read(unit,'(8I6)') pcell%symrel
299   read(unit,'(a)')
300   read(unit,'(8I6)') pcell%typat
301 
302 ! read in real
303   read(unit,'(a)')
304   read(unit,'(3E20.10)') pcell%acell
305   read(unit,'(a)')
306   read(unit,'(3E20.10)') pcell%amu
307   read(unit,'(a)')
308   read(unit,'(3E20.10)') pcell%dielt
309   read(unit,'(a)')
310   read(unit,'(3E20.10)') pcell%dyewq0
311   read(unit,'(a)')
312   read(unit,'(3E20.10)') pcell%gmet
313   read(unit,'(a)')
314   read(unit,'(3E20.10)') pcell%gprim
315   read(unit,'(a)')
316   read(unit,'(3E20.10)') pcell%rcan
317   read(unit,'(a)')
318   read(unit,'(3E20.10)') pcell%rmet
319   read(unit,'(a)')
320   read(unit,'(3E20.10)') pcell%rprim
321   read(unit,'(a)')
322   read(unit,'(3E20.10)') pcell%rpt
323   read(unit,'(a)')
324   read(unit,'(3E20.10)') pcell%trans
325   read(unit,'(a)')
326   read(unit,'(3E20.10)') pcell%wghatm
327   read(unit,'(a)')
328   read(unit,'(3E20.10)') pcell%xred
329   read(unit,'(a)')
330   read(unit,'(3E20.10)') pcell%zeff
331 
332   close(unit)
333   
334  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-2018 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

PARENTS

CHILDREN

SOURCE

367 subroutine write_primcell_ddb_info (filename,pcell)
368  
369  use defs_basis
370 
371 !This section has been created automatically by the script Abilint (TD).
372 !Do not modify the following lines by hand.
373 #undef ABI_FUNC
374 #define ABI_FUNC 'write_primcell_ddb_info'
375 !End of the abilint section
376 
377  implicit none
378 
379 !Arguments ------------------------------------
380 
381  character(len=*), intent(in) :: filename
382  type(primcell_ddb_info), intent(in) :: pcell
383 
384 !Local variables-------------------------------
385  integer :: unit
386  character(len=500) :: msg
387 
388 ! *************************************************************************
389 
390   if (open_file(filename,msg,newunit=unit, form="formatted",status="unknown") /= 0) then
391     MSG_ERROR(msg)
392   end if
393 
394 ! read out dimensions
395   write(unit,'(a,I6)') 'pcell%brav   ', pcell%brav
396   write(unit,'(a,I6)') 'pcell%dipdip ', pcell%dipdip
397   write(unit,'(a,I6)') 'pcell%mpert  ', pcell%mpert
398   write(unit,'(a,I6)') 'pcell%msym   ', pcell%msym
399   write(unit,'(a,I6)') 'pcell%natom  ', pcell%natom
400   write(unit,'(a,I6)') 'pcell%nrpt   ', pcell%nrpt
401   write(unit,'(a,I6)') 'pcell%nsym   ', pcell%nsym
402   write(unit,'(a,I6)') 'pcell%ntypat ', pcell%ntypat
403 
404 ! write out scalar reals
405   write(unit,'(a,E20.10)') 'pcell%ucvol  ', pcell%ucvol
406 
407 ! write out int
408   write(unit,'(a)') 'pcell%indsym'
409   write(unit,'(8I6)') pcell%indsym
410   write(unit,'(a)') 'pcell%symrec'
411   write(unit,'(8I6)') pcell%symrec
412   write(unit,'(a)') 'pcell%symrel'
413   write(unit,'(8I6)') pcell%symrel
414   write(unit,'(a)') 'pcell%typat'
415   write(unit,'(8I6)') pcell%typat
416 
417 ! write out real
418   write(unit,'(a)') 'pcell%acell'
419   write(unit,'(3E20.10)') pcell%acell
420   write(unit,'(a)') 'pcell%amu'
421   write(unit,'(3E20.10)') pcell%amu
422   write(unit,'(a)') 'pcell%dielt'
423   write(unit,'(3E20.10)') pcell%dielt
424   write(unit,'(a)') 'pcell%dyewq0'
425   write(unit,'(3E20.10)') pcell%dyewq0
426   write(unit,'(a)') 'pcell%gmet'
427   write(unit,'(3E20.10)') pcell%gmet
428   write(unit,'(a)') 'pcell%gprim'
429   write(unit,'(3E20.10)') pcell%gprim
430   write(unit,'(a)') 'pcell%rcan'
431   write(unit,'(3E20.10)') pcell%rcan
432   write(unit,'(a)') 'pcell%rmet'
433   write(unit,'(3E20.10)') pcell%rmet
434   write(unit,'(a)') 'pcell%rprim'
435   write(unit,'(3E20.10)') pcell%rprim
436   write(unit,'(a)') 'pcell%rpt'
437   write(unit,'(3E20.10)') pcell%rpt
438   write(unit,'(a)') 'pcell%trans'
439   write(unit,'(3E20.10)') pcell%trans
440   write(unit,'(a)') 'pcell%wghatm'
441   write(unit,'(3E20.10)') pcell%wghatm
442   write(unit,'(a)') 'pcell%xred'
443   write(unit,'(3E20.10)') pcell%xred
444   write(unit,'(a)') 'pcell%zeff'
445   write(unit,'(3E20.10)') pcell%zeff
446 
447   close(unit)
448   
449  end subroutine write_primcell_ddb_info