TABLE OF CONTENTS


ABINIT/mksupercell [ Functions ]

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

COPYRIGHT

  Copyright (C) 1998-2018 ABINIT group (DJA)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .
  For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

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

 37 #if defined HAVE_CONFIG_H
 38 #include "config.h"
 39 #endif
 40 
 41 #include "abi_common.h"
 42 
 43 subroutine mksupercell(xred_org,magv_org,rprimd_org,nat_org,nat_sc,xred_sc,magv_sc,rprimd_sc,ext,prtvol) 
 44 
 45  use defs_basis
 46  use m_profiling_abi
 47 
 48 !This section has been created automatically by the script Abilint (TD).
 49 !Do not modify the following lines by hand.
 50 #undef ABI_FUNC
 51 #define ABI_FUNC 'mksupercell'
 52 !End of the abilint section
 53 
 54  implicit none
 55 
 56 !Arguments ------------------------------------
 57 !scalars
 58  integer,intent(in)              :: nat_org,nat_sc
 59  integer,intent(in),optional     :: prtvol
 60 !arrays
 61  real(dp),intent(in)             :: rprimd_org(3,3)
 62  integer,intent(in)              :: ext(3)
 63  real(dp),intent(in)             :: xred_org(3,nat_org)
 64  real(dp),intent(out)            :: xred_sc(3,nat_sc) 
 65  real(dp),intent(out)            :: magv_sc(nat_sc)
 66  real(dp),intent(out)            :: rprimd_sc(3,3)
 67  integer,intent(in),optional     :: magv_org(nat_org)
 68 
 69 
 70 !Local variables-------------------------------
 71 !scalars
 72  integer                      :: prtvoll,ix,iy,iz,nprcl,iprcl,jdim,iatom
 73 !arrays
 74  real(dp)                     :: magvv_org(nat_org)
 75  real(dp),allocatable         :: transv(:,:,:)
 76 
 77 ! *************************************************************************
 78 
 79 !DEBUG
 80 !write(std_out,*)'mksupercell: enter'
 81 !END DEBUG
 82 
 83  if (present(magv_org)) then
 84    magvv_org=magv_org
 85  else
 86    magvv_org=(/ (1, iatom=1,nat_org)  /) 
 87  end if
 88 
 89  if (present(prtvol)) then
 90    prtvoll=prtvol
 91  else
 92    prtvoll=1
 93  end if
 94 
 95  rprimd_sc=reshape((/ (rprimd_org(ix,:)*ext(ix) ,ix=1,3) /),(/3,3 /))
 96  nprcl=product(ext)
 97  ABI_ALLOCATE(transv,(3,nat_org,nprcl))
 98 
 99  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/) )
100 
101 !DEBUG
102 !write(std_out,*)'mksupercell: xred_org ' ,xred_org
103 !END DEBUG
104 
105  do iprcl=1,nprcl
106    xred_sc(:,1+(iprcl-1)*nat_org:iprcl*nat_org)=xred_org+transv(:,:,iprcl)
107    magv_sc(1+(iprcl-1)*nat_org:iprcl*nat_org)=magv_org
108  end do
109 
110 
111  do jdim=1,3
112    xred_sc(jdim,:)=xred_sc(jdim,:)/ext(jdim)
113  end do
114 
115 !DEBUG
116 !write(std_out,*)'mksupercell: xred_sc ', xred_sc
117 !write(std_out,*)'mksupercell: magv_sc ', magv_sc
118 !END DEBUG
119 
120  ABI_DEALLOCATE(transv)
121  
122 !DEBUG
123 !write(std_out,*)'mksupercell: leave'
124 !END DEBUG
125 end subroutine mksupercell