TABLE OF CONTENTS


ABINIT/m_supercell [ Modules ]

[ Top ] [ Modules ]

NAME

 m_supercell

FUNCTION

 Module for using a supercell, in particular for phonon displacement freezing.
 Container type is defined, and destruction, print subroutines
 as well as the central init_supercell

COPYRIGHT

 Copyright (C) 2010-2018 ABINIT group (MJV, DJA)
 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 module m_supercell
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32 
33  use m_symtk,    only : matr3inv
34  use m_copy,     only : alloc_copy
35  use m_io_tools, only : open_file
36  use m_fstrings, only : int2char4, write_num
37 
38  implicit none
39 
40  private

defs_abitypes/supercell_type [ Types ]

[ Top ] [ defs_abitypes ] [ Types ]

NAME

 supercell_type

FUNCTION

 structure for a supercell constructed from a basic rprimd and xcart, with indexing to original atoms
 The supercell may not be oriented the same way as the original cell, if you can reduce it by symmetry

SOURCE

54  type, public :: supercell_type
55    integer :: natom_primcell                        ! number of atoms in primitive cell
56    integer :: natom                                 ! number of atoms in supercell
57    integer :: ntypat                                ! number of atom types
58    integer :: ncells                                ! number of unit cells in supercell
59    integer :: rlatt(3,3)                            ! matrix for multiplicity of supercell
60    real(dp) :: rprimd(3,3)                          ! new lattice vectors for supercell
61    real(dp) :: qphon(3)                             ! phonon q vector used to generate scell, if any
62    real(dp), allocatable :: xcart(:,:)              ! (3, natom) positions of atoms
63    real(dp), allocatable :: xcart_ref(:,:)          ! (3, natom) equilibrium positions of atoms
64    integer, allocatable :: atom_indexing(:)         ! (natom) indexes original atom: 1..natom_primcell
65    integer, allocatable :: uc_indexing(:,:)         ! (3, natom) indexes unit cell atom is in:
66    integer, allocatable :: typat(:)                 ! (3, natom) positions of atoms
67    real(dp), allocatable :: znucl(:)                ! (ntypat) nuclear charges of species
68    integer, allocatable :: rvecs(:,:)               ! supercell vectors
69 
70  end type supercell_type
71 
72  public :: init_supercell_for_qpt
73  public :: init_supercell
74  public :: freeze_displ_supercell
75  public :: prt_supercell_for_qpt
76  public :: prt_supercell
77  public :: copy_supercell
78  public :: distance_supercell
79  public :: findBound_supercell
80  public :: getPBCIndexes_supercell
81  public :: destroy_supercell
82 
83  public :: mksupercell  !  computes atomic positons, magnetic ordering of supercell

m_effective_potential/getPBCIndexes_supercell [ Functions ]

[ Top ] [ m_effective_potential ] [ Functions ]

NAME

FUNCTION

 Get the index of the cell by using PBC

INPUTS

 index  = index of the cell into the supercell
 ncell = number of total cell

OUTPUT

 index  = index of the cell into the supercell with PBC

SOURCE

697 subroutine getPBCIndexes_supercell(index,ncell)
698 
699 
700 !This section has been created automatically by the script Abilint (TD).
701 !Do not modify the following lines by hand.
702 #undef ABI_FUNC
703 #define ABI_FUNC 'getPBCIndexes_supercell'
704 !End of the abilint section
705 
706  implicit none
707 
708 !Arguments ---------------------------------------------
709   integer, intent(inout)  :: index(3)
710   integer, intent(in) :: ncell(3)
711 !Local variables ---------------------------------------
712   integer :: ii
713 ! *********************************************************************
714 
715   do ii=1,3
716     do while (index(ii) > ncell(ii))
717       index(ii) = index(ii) - ncell(ii)
718     end do
719     do while (index(ii) <= 0)
720       index(ii) = index(ii) + ncell(ii)
721     end do
722   end do
723 
724 end subroutine getPBCIndexes_supercell

m_supercell/distance_supercell [ Functions ]

[ Top ] [ m_supercell ] [ Functions ]

NAME

FUNCTION

 compute the distance_supercell betwen 2 atoms in different cell

INPUTS

 xcart1(3) = cartesian coordinates of the first atom
 xcart1(3) = cartesian coordinates of the second atom
 rprimd(3,3) = primitive lattice vectors
 cell1(3) = index of the cell of the first atom (for example -1 0 2)
 cell2(3) = index of the cell of the second atom (for example  0 0 2)

OUTPUT

 distance_supercell = distance_supercell between the 2 atoms

SOURCE

794 function distance_supercell(xcart1,xcart2,rprimd,cell1,cell2) result(dist)
795 
796 !Arguments ------------------------------------
797 !scalar
798 !array
799 
800 !This section has been created automatically by the script Abilint (TD).
801 !Do not modify the following lines by hand.
802 #undef ABI_FUNC
803 #define ABI_FUNC 'distance_supercell'
804 !End of the abilint section
805 
806   real(dp),intent(in):: rprimd(3,3)
807   real(dp),intent(in):: xcart1(3),xcart2(3)
808   integer,intent(in) :: cell1(3),cell2(3)
809   real(dp) :: dist
810 !Local variables -------------------------------
811   real(dp) :: rpt1(3),rpt2(3)
812   integer  :: mu
813 !! *************************************************************************
814   do mu=1,3
815     rpt1(mu) = cell1(1)*rprimd(mu,1)+cell1(2)*rprimd(mu,2)+cell1(3)*rprimd(mu,3)
816     rpt2(mu) = cell2(1)*rprimd(mu,1)+cell2(2)*rprimd(mu,2)+cell2(3)*rprimd(mu,3)
817   end do
818 
819   dist = ((xcart2(1)+rpt2(1)-xcart1(1)-rpt1(1))**2+&
820 &         (xcart2(2)+rpt2(2)-xcart1(2)-rpt1(2))**2+&
821 &         (xcart2(3)+rpt2(3)-xcart1(3)-rpt1(3))**2)**0.5
822 
823 end function distance_supercell

m_supercell/freeze_displ_supercell [ Functions ]

[ Top ] [ m_supercell ] [ Functions ]

NAME

 freeze_displ_supercell

FUNCTION

 Freeze a specific displacement phonon field into the supercell scell

INPUTS

 displ = phonon displacement vectors for this mode
 freeze_displ = desired amplitude for phonon displacement along displ.
    for thermal displacement use sqrt[ (1/2 + bose_einstein(freq,T)) / freq ]
 scell = supercell structure with reference atomic positions etc...

OUTPUT

 scell = supercell structure: xcart will be updated with phonon displacement

PARENTS

      freeze_displ_allmodes,m_phonons

CHILDREN

SOURCE

389 subroutine freeze_displ_supercell (displ,freeze_displ,scell)
390 
391 
392 !This section has been created automatically by the script Abilint (TD).
393 !Do not modify the following lines by hand.
394 #undef ABI_FUNC
395 #define ABI_FUNC 'freeze_displ_supercell'
396 !End of the abilint section
397 
398  implicit none
399 
400 !Arguments ------------------------------------
401 !scalars
402  type(supercell_type), intent(inout) :: scell
403  real(dp), intent(in) :: freeze_displ
404 !arrays
405  real(dp), intent(in) :: displ(2,3*scell%natom_primcell)
406 
407 !Local variables-------------------------------
408 !scalar
409  integer :: iatom, ipratom
410  complex(dpc) :: expqdotr, j=cmplx(zero,one)
411  complex(dpc) :: phase
412  complex(dpc) :: zdispl(3,scell%natom_primcell)
413 ! *************************************************************************
414 
415  zdispl = (cmplx(reshape(displ(1,:), (/3,scell%natom_primcell/)),&
416 &                reshape(displ(2,:), (/3,scell%natom_primcell/))))
417 
418  ! fix gauge by imposing real displacement for first atom in first direction
419  ! multiply by normalized complex conjugate of first element
420  ! NB 6 March 2018: this may be imposing a positive (not just real) displacement for 1st atom along x!!! 
421  ! That might be problematic below, though for the thermal displacement method freeze_displ swaps sign for each new mode 
422  phase = cmplx(one,zero)
423  if (abs(zdispl(1,1)) > tol10) then
424    phase = conjg(zdispl(1,1)) / abs(zdispl(1,1))
425  end if
426 
427  do iatom = 1, scell%natom
428    expqdotr = exp(j*two_pi*(scell%qphon(1)*scell%uc_indexing(1,iatom) &
429 &                          +scell%qphon(2)*scell%uc_indexing(2,iatom) &
430 &                          +scell%qphon(3)*scell%uc_indexing(3,iatom)))
431 
432 ! this is offset in zdispl vector due to primitive cell atom position
433    ipratom = scell%atom_indexing(iatom)
434 
435 !add real part of displacement times Bloch phase
436    scell%xcart(:,iatom) = scell%xcart(:,iatom) &
437 &        + freeze_displ * real(expqdotr * zdispl(:,ipratom) * phase)
438 
439 !   scell%xcart(:,iatom) = scell%xcart(:,iatom) &
440 !&        + freeze_displ * cos(qdotr) * displ(1,ipratom+1:ipratom+3) &
441 !&        - freeze_displ * sin(qdotr) * displ(2,ipratom+1:ipratom+3)
442  end do
443 
444 end subroutine freeze_displ_supercell

m_supercell/init_supercell [ Functions ]

[ Top ] [ m_supercell ] [ Functions ]

NAME

 init_supercell

FUNCTION

 Initialize scell structure, from unit cell vectors, and atoms, based on rlatt multiplicity matrix

INPUTS

 natom_primcell = number of atoms in primitive cell
 rlatt(3,3) = multiplicity of primtive unit cells in supercell
 rprimd_primcell(3,3) = real space lattice vectors (bohr)
 typat_primcell(natom) = types of all atoms in primitive cell
 xcart_primcell(3,natom) = cartesian positions of atoms in primitive cell
 znucl = nuclear charges for all species

 ordering = if true,  typat will be 1 1 1 1 1 1 2 2 2 2 2 2 3 3 3 3 ....
            if false, typat will be 1 2 3 4 1 2 3 4 1 2 3 4 1 2 3 4 ....

OUTPUT

 scell = supercell structure to be initialized

PARENTS

      m_effective_potential,m_fit_polynomial_coeff,m_phonons,m_supercell

CHILDREN

SOURCE

215 subroutine init_supercell(natom_primcell, rlatt, rprimd_primcell, typat_primcell, xcart_primcell, znucl, scell, ordering)
216 
217 
218 !This section has been created automatically by the script Abilint (TD).
219 !Do not modify the following lines by hand.
220 #undef ABI_FUNC
221 #define ABI_FUNC 'init_supercell'
222 !End of the abilint section
223 
224  implicit none
225 !Arguments ------------------------------------
226 !scalars
227  integer, intent(in) :: natom_primcell
228  type(supercell_type), intent(out) :: scell
229  logical,optional,intent(in) :: ordering
230 !arrays
231  integer , intent(in) :: rlatt(3,3)
232  integer , intent(in) :: typat_primcell(natom_primcell)
233  real(dp), intent(in) :: znucl(:)
234  real(dp), intent(in) :: rprimd_primcell(3,3)
235  real(dp), intent(in) :: xcart_primcell(3,natom_primcell)
236 
237 !local
238 !scalars
239  integer :: iatom_supercell, i1,i2,i3, iatom, icell
240 
241 !arrays
242 
243  scell%natom_primcell = natom_primcell
244  scell%rlatt = rlatt
245  scell%ncells = rlatt(1,1)*rlatt(2,2)*rlatt(3,3)
246  scell%rprimd(:,1) = rprimd_primcell(:,1) * rlatt(1,1)
247  scell%rprimd(:,2) = rprimd_primcell(:,2) * rlatt(2,2)
248  scell%rprimd(:,3) = rprimd_primcell(:,3) * rlatt(3,3)
249 
250  scell%ntypat = size(znucl)
251  ABI_ALLOCATE(scell%znucl,(scell%ntypat))
252  scell%znucl(:) = znucl(:)
253 
254 !number of atoms in full supercell
255  scell%natom= natom_primcell*scell%ncells
256  ABI_ALLOCATE(scell%xcart,(3,scell%natom))
257  ABI_ALLOCATE(scell%xcart_ref,(3,scell%natom))
258  ABI_ALLOCATE(scell%typat,(scell%natom))
259  ABI_ALLOCATE(scell%atom_indexing,(scell%natom))
260  ABI_ALLOCATE(scell%uc_indexing,(3,scell%natom))
261  ABI_ALLOCATE(scell%rvecs, (3, scell%ncells))
262 
263  iatom_supercell = 0
264  icell =0
265  do i1 = 1, rlatt(1,1)
266    do i2 = 1, rlatt(2,2)
267      do i3 = 1, rlatt(3,3)
268 
269        icell=icell+1
270        scell%rvecs(:,icell)=(/i1-1, i2-1, i3-1/)
271 
272 
273        do iatom = 1, natom_primcell
274          iatom_supercell = iatom_supercell + 1
275          scell%uc_indexing(:,iatom_supercell) = (/i1-1,i2-1,i3-1/)
276          scell%xcart_ref(:,iatom_supercell) = xcart_primcell(:,iatom) &
277 &            + matmul(rprimd_primcell,scell%uc_indexing(:,iatom_supercell))
278          scell%atom_indexing(iatom_supercell) = iatom
279          scell%typat(iatom_supercell) = typat_primcell(iatom)
280        end do
281      end do
282    end do
283  end do
284 
285  ABI_CHECK(iatom_supercell == scell%natom, "iatom_supercell /= scell%natom")
286  if(iatom_supercell /= scell%natom) then
287     print *, "iatom_supercell /= scell%natom"
288  endif
289 
290 
291 
292  scell%xcart = scell%xcart_ref
293  scell%qphon = zero
294 
295  if (present(ordering)) then
296    if (ordering) then
297      call order_supercell_typat(scell)
298    end if
299  end if
300 
301 end subroutine init_supercell

m_supercell/init_supercell_for_qpt [ Functions ]

[ Top ] [ m_supercell ] [ Functions ]

NAME

 init_supercell_for_qpt

FUNCTION

 Initialize scell structure, from unit cell vectors, and atoms, based on qpoint chosen

INPUTS

 natom_primcell = number of atoms in primitive cell
 qphon(3) = phonon wavevector
      find smallest supercell which will accomodate phonon qphon = (1/2,1/2,1/2)
 rprimd_primcell(3,3) = real space lattice vectors (bohr)
 typat_primcell = types of atoms
 xcart_primcell(3,natom) = cartesian positions of atoms in primitive cell
 znucl = nuclear charges for all species
 ordering = if true,  typat will be 1 1 1 1 1 1 2 2 2 2 2 2 3 3 3 3 ....
            if false, typat will be 1 2 3 4 1 2 3 4 1 2 3 4 1 2 3 4 ....

OUTPUT

 scell = supercell structure to be initialized

PARENTS

      freeze_displ_allmodes

CHILDREN

SOURCE

117 subroutine init_supercell_for_qpt(natom_primcell, qphon, rprimd_primcell, &
118 &    typat_primcell, xcart_primcell, znucl, scell, ordering)
119 
120 
121 !This section has been created automatically by the script Abilint (TD).
122 !Do not modify the following lines by hand.
123 #undef ABI_FUNC
124 #define ABI_FUNC 'init_supercell_for_qpt'
125 !End of the abilint section
126 
127  implicit none
128 !Arguments ------------------------------------
129 !scalars
130  integer, intent(in) :: natom_primcell
131  type(supercell_type), intent(out) :: scell
132  logical,optional,intent(in) :: ordering
133 !arrays
134  integer , intent(in) :: typat_primcell(natom_primcell)
135  real(dp), intent(in) :: qphon(3)
136  real(dp), intent(in) :: znucl(:)
137  real(dp), intent(in) :: rprimd_primcell(3,3)
138  real(dp), intent(in) :: xcart_primcell(3,natom_primcell)
139 
140 !Local variables-------------------------------
141 !scalar
142  integer :: ii, maxsc, iscmult
143  real(dp) :: qbymult
144 !arrays
145  integer :: rlatt(3,3) ! number of primitive cells in each direction for the supercell
146  character(len=500) :: msg
147 ! *************************************************************************
148 
149 
150 ! maximum number of unit cells in a given direction
151  maxsc = 10
152 
153 ! find smallest supercell which will accomodate phonon.
154 ! FIXME: for the moment, just get smallest multiple along each direction, with an upper bound
155  rlatt = 0
156  rlatt(1,1) = -1
157  rlatt(2,2) = -1
158  rlatt(3,3) = -1
159  do ii=1,3
160    do iscmult=1,maxsc
161      qbymult = qphon(ii)*iscmult
162      if (abs(qbymult - int(qbymult)) < tol10) then
163        rlatt(ii,ii) = iscmult
164        exit
165      end if
166    end do
167    if (rlatt(ii,ii) == -1) then
168      write(msg,'(a,I4,a,I7,2a,3E20.10)')' No supercell found with less than ', &
169 &           maxsc,' unit cells in direction ', &
170 &           ii, ch10, ' qphon = ', qphon
171      MSG_ERROR(msg)
172    end if
173  end do
174 
175  if (present(ordering)) then
176    call init_supercell(natom_primcell, rlatt, rprimd_primcell, typat_primcell, xcart_primcell, znucl, scell, ordering)
177  else
178    call init_supercell(natom_primcell, rlatt, rprimd_primcell, typat_primcell, xcart_primcell, znucl, scell)
179  end if
180 
181  scell%qphon = qphon
182 
183 end subroutine init_supercell_for_qpt

m_supercell/mksupercell [ Functions ]

[ Top ] [ m_supercell ] [ Functions ]

NAME

  mksupercell

FUNCTION

  computes atomic positons, magnetic ordering of supercell

INPUTS

  magv_org (optional) magnetic ordering of atoms in primitive cell,
   ordering of atoms given als 1 and -1, if not given fm is assumed
  xred_org relative position of atoms in primitive cell
  rprimd_org unit cell dimensions of primitive cell
  natom=number of atoms in unit cell
  option= 1 output ion-ion distances / 2 output ordering of ion-ion distances / 3 output variables in varlist
           according to ion-ion distances * magnetic ordering

OUTPUT

  magv_sc magnetic ordering of atoms in supercell
  xred_sc relative position of atoms in supercell
  rprimd_sc unit cell dimensions of supercell

PARENTS

      pawuj_det

CHILDREN

SOURCE

918 subroutine mksupercell(xred_org,magv_org,rprimd_org,nat_org,nat_sc,xred_sc,magv_sc,rprimd_sc,ext,prtvol)
919 
920 
921 !This section has been created automatically by the script Abilint (TD).
922 !Do not modify the following lines by hand.
923 #undef ABI_FUNC
924 #define ABI_FUNC 'mksupercell'
925 !End of the abilint section
926 
927  implicit none
928 
929 !Arguments ------------------------------------
930 !scalars
931  integer,intent(in)              :: nat_org,nat_sc
932  integer,intent(in),optional     :: prtvol
933 !arrays
934  real(dp),intent(in)             :: rprimd_org(3,3)
935  integer,intent(in)              :: ext(3)
936  real(dp),intent(in)             :: xred_org(3,nat_org)
937  real(dp),intent(out)            :: xred_sc(3,nat_sc)
938  real(dp),intent(out)            :: magv_sc(nat_sc)
939  real(dp),intent(out)            :: rprimd_sc(3,3)
940  integer,intent(in),optional     :: magv_org(nat_org)
941 
942 !Local variables-------------------------------
943 !scalars
944  integer :: prtvoll,ix,iy,iz,nprcl,iprcl,jdim,iatom
945 !arrays
946  real(dp) :: magvv_org(nat_org)
947  real(dp),allocatable :: transv(:,:,:)
948 
949 ! *************************************************************************
950 
951  if (present(magv_org)) then
952    magvv_org=magv_org
953  else
954    magvv_org=(/ (1, iatom=1,nat_org)  /)
955  end if
956 
957  if (present(prtvol)) then
958    prtvoll=prtvol
959  else
960    prtvoll=1
961  end if
962 
963  rprimd_sc=reshape((/ (rprimd_org(ix,:)*ext(ix) ,ix=1,3) /),(/3,3 /))
964  nprcl=product(ext)
965  ABI_ALLOCATE(transv,(3,nat_org,nprcl))
966 
967  transv=reshape((/ (((((/ ix,iy,iz /),iatom=1,nat_org),ix=0,ext(1)-1),iy=0,ext(2)-1),iz=0,ext(3)-1) /), (/ 3, nat_org,nprcl/) )
968 
969 !DEBUG
970 !write(std_out,*)'mksupercell: xred_org ' ,xred_org
971 !END DEBUG
972 
973  do iprcl=1,nprcl
974    xred_sc(:,1+(iprcl-1)*nat_org:iprcl*nat_org)=xred_org+transv(:,:,iprcl)
975    magv_sc(1+(iprcl-1)*nat_org:iprcl*nat_org)=magv_org
976  end do
977 
978 
979  do jdim=1,3
980    xred_sc(jdim,:)=xred_sc(jdim,:)/ext(jdim)
981  end do
982 
983 !DEBUG
984 !write(std_out,*)'mksupercell: xred_sc ', xred_sc
985 !write(std_out,*)'mksupercell: magv_sc ', magv_sc
986 !END DEBUG
987 
988  ABI_DEALLOCATE(transv)
989 
990 end subroutine mksupercell

m_supercell/order_supercell_typat [ Functions ]

[ Top ] [ m_supercell ] [ Functions ]

NAME

 order_supercell_typat

FUNCTION

 Re-order atoms in place for types

INPUTS

 scell = supercell structure with reference atomic positions etc...

OUTPUT

 scell = supercell structure: typat, xcart and so on will be updated

PARENTS

      m_supercell

CHILDREN

SOURCE

325 subroutine order_supercell_typat (scell)
326 
327 
328 !This section has been created automatically by the script Abilint (TD).
329 !Do not modify the following lines by hand.
330 #undef ABI_FUNC
331 #define ABI_FUNC 'order_supercell_typat'
332 !End of the abilint section
333 
334  implicit none
335 
336 !Arguments ------------------------------------
337 !scalars
338  type(supercell_type), intent(inout) :: scell
339 
340 ! local tmp variables
341  integer :: itypat, iatom_supercell, iatom
342  type(supercell_type) :: scell_tmp
343 
344  call copy_supercell (scell,scell_tmp)
345 
346  iatom_supercell = 0
347  do itypat = 1, scell%ntypat
348    do iatom = 1, scell%natom
349      if (scell_tmp%typat(iatom) /= itypat) cycle
350      iatom_supercell = iatom_supercell + 1
351      scell%xcart(:,iatom_supercell) = scell_tmp%xcart(:,iatom)
352      scell%xcart_ref(:,iatom_supercell) = scell_tmp%xcart_ref(:,iatom)
353      scell%atom_indexing(iatom_supercell) = scell_tmp%atom_indexing(iatom)
354      scell%uc_indexing(:,iatom_supercell) = scell_tmp%uc_indexing(:,iatom)
355      scell%typat(iatom_supercell) = scell_tmp%typat(iatom)
356    end do
357  end do
358 
359  call destroy_supercell(scell_tmp)
360 
361 end subroutine order_supercell_typat